Feature reconstruction based on t-SNE: an approach for fault diagnosis of rotating machinery

2.1. LCD

The LCD method can self-adaptively decompose a complex signal into a sum of ISCs and a residue, with the assumption that any two ISCs are independent of each other. The decomposition process is shown in Fig. 2. In the figure, the standard deviation (SD) is defined as the criterion to determine whether the decomposition ends. In this study, the value of the SD is set as 0.3. Therefore, when the SD is smaller than 0.3, the decomposition process ends.

The ISCs must meet the following two conditions:

1) In the entire data set, the signal is monotonic between any two adjacent extreme points;

2) Among the whole data set, suppose the extrema of the signal $x\left(t\right)$ are ($t_k$, $x_k$) ($k=$ 1, 2, …, $M$). Thus, any two adjacent maxima (minima), ($t_{k-1}$, $x_{k-1}$) and ($t_{k+1}$, $x_{k+1}$), are connected by a straight line. For the intermediate minima (maxima) $(t_k,x_k)$, its corresponding point ($t_k$, $A_k$) on this straight line is as follows:

To ensure the smoothness and symmetry of the ISC, the proportions of $x_k$ and $A_k$remain constant, as described below:

Fig. 2Decomposition process of the LCD

Generally, $a$ is set as 0.5; thus, $A_k/x_k=$ –1. At the time, $x_k$ and $A_k$ are symmetrical about the $x$-axis, thus ensuring the symmetry and smoothness of the ISC component waveform [13].

Based on the definition of the ISC component, a complex signal $x\left(t\right)$ ($t>$ 0) can be decomposed into a sum of ISCs using the LCD method as follows:

1) Suppose the number of extrema of the signal $x\left(t\right)$ is $M$, and determine all the extrema ($t_k$, $x_k$) ($k=$ 1, …, $M$) of the signal $x\left(t\right)$.

2) Compute $A_k$($k=$ 1, …, $M$) according to Eq. (1). Based on the following equation, the corresponding $L_k$($k=$ 1, …, $M$) can be calculated:

To obtain the values of $L_1$ and $L_M$, we should extend the boundary of the data; this process can be completed via different methods [14-16]. By extension, two end extrema can be obtained and written as ($t_0$, $x_0$) and ($t_{M+1}$, $x_{M+1}$) [17]. Using Eq. (1) and Eq. (3), we can obtain $L_1$ and $L_M$.

3) All the $L_k$($k=$ 1, …, $M$) are connected by a cubic spline line $S{L}_1\left(t\right)$, which is defined as the baseline of the LCD.

4) The difference between the data and the baseline $S{L}_1\left(t\right)$ is the first component $h_1\left(t\right)$:

5) Separate the first $IS{C}_1$ from the initial data$x\left(t\right)$:

6) Take the residue $r_1$ as the original signal to be processed. Repeat the above process until the residue is constant, a monotonic function, or a function with no more than twenty extrema.

$IS{C}_1$, …, $IS{C}_n$ and the residue are obtained as follows:

where $c_i\left(t\right)$ is the $i$th ISC, and $r_n\left(t\right)$ is the residue.

To illustrate the validity of the LCD method, a sample vibration signal of a bearing and its ISCs decomposed by LCD are shown in Fig. 3. The figure shows that the original vibration signal can be effectively decomposed into 9 ISCs.

Fig. 3An original vibration signal of a bearing and its ICSs decomposed by LCD

a)

b)

c)

d)

e)

f)

g)

h)

i)

j)

2.2. SampEn

The LCD method is applied to decompose the original signal into several ISCs with different time scales and characteristics. Next, the SampEn is used to characterize the probability of the new model of the signal, i.e., the change in the values of SampEn when the centrifugal pump works under different conditions can reflect the different fault modes of the centrifugal pump. Thus, the LCD-SampEn can be applied to obtain the feature vector.

The SampEn of a signal of length $N$($x_1$, $x_2$, …, $x_N$) is defined as:

where:

And:

In Eq. (8), $N$ denotes the number of data points, $m$ specifies the pattern length, and $\tau$ is the time delay. In Eq. (4), the symbol $\mathrm{\Theta }$ is the Heaviside function, given as:

The distance $d_{ij}^m$ between $X_i^m$ and $X_j^m$ ($m$-dimensional pattern vector) is defined as:

And:

The symbol $r$ represents a predetermined tolerance value, which is defined as:

where $k$ is a constant ($k>$ 0), and $std$($\cdot$) represents the standard deviation of the signal. The two patterns $i$ and $j$ of $m$ measurements of the signal are similar if the distance $d_{ij}^m$ between any pair of corresponding measurements of $X_i^m$ is less than or equal to $r$ [18].

2.3. t-SNE

Stochastic neighbor embedding (SNE) is one of the best performing nonlinear manifold learning algorithms and whose core idea is to maintain the consistency of the neighborhood probability distribution between high-dimensional and low-dimensional space. SNE transfers the traditional Euclidean distance-based similarity measurement to the conditional probability-based similarity measurement: in high-dimensional observation space, the Gaussian distribution is adopted to simulate the similarity relationship between observation samples. The similarity between $x_i$ and $x_j$$p_{j|i}$ is denoted as follows:

where ${\sigma }_i$ is the bandwidth of the Gaussian kernel function in the observation sample $x_i$. $p_{j|i}$ is the probability of $x_j$ chosen by $x_i$ as its neighbor sample. The parameter $p_{j|i}$ obeys a Gaussian distribution, for which the variance is ${\sigma }_i^2$ and the mean value is $x_i$. The probability of $x_i$ and $x_j$ being adjacent to each other is:

In low-dimensional space, SNE continues to adopt the Gaussian distribution to measure the similarity between low-distinction samples. However, two obvious shortcomings exist with SNE: (1) the objective function is too complex to optimize, and the gradient is not as concise as desired; (2) the so-called “crowding problem”; that is, when the data are far apart from each other in high-dimensional space, they must be gathered in the process of mapping to low-dimensional space. In response, t-distributed stochastic neighbor embedding (t-SNE) is introduced to alleviate these problems [19].

To solve the first problem, the characteristic of symmetry is adopted to simplify the objective function and optimize the gradient form, which is referred to as symmetric SNE. According to probability theory, the SNE objective function minimizes the sum of the distances of the conditional probability distribution $p_{j|i}$ (high-dimension) and $p_{i|j}$ (low-dimension) for the corresponding points. It is equal to the following joint probability distributions of $P$ (high-dimension) and $Q$ (low-dimension):

After adopting a joint probability distribution instead of a conditional probability distribution, the formula is more concise and understandable.

To solve the second problem, the t-distribution function is introduced to alleviate the “crowding problem”. In other words, the t-distribution function is used to measure the similarity of points in low-dimensional space. The joint probability distribution function is as follows:

Here, the t-distribution function (DOF is 1) is applied because of its special advantageous characteristic: $(1+{‖y_i-y_j‖}_2^2{)}^{-1}$ is the reciprocal of the distance of points far from each other in low-dimensional space to ${‖y_i-y_k‖}_{}^2$. In other words, in low-dimensional space, the presentation of the joint probability distribution is insensitive to the distance of points. In addition, in theory, the t-distribution function offers the same performance as the Gaussian function because the t-distribution function can be expressed as the infinity Gaussian function. Thus, the gradient of t-SNE is:

In conclusion, t-SNE focuses on the following:

1) The characteristic of non-similarity is associated with points far from each other. The characteristic of similarity is associated with points close to each other. The t-distribution function is introduced to exert an “exclusion” to process the non-similar points;

2) The application of the t-distribution function makes it easier to optimize.

The t-SNE approach applies a probability density distribution function to measure the similarity and distribution characteristics. Next, t-SNE accomplishes property preservation by minimizing the K-L distance probability density in both high-dimensional and low-dimensional space. The t-distribution function is employed to measure the similarity of points in low-dimensional space, thereby simplifying the gradient form and improves computational speed. Most importantly, the improved processing of the non-similarity points alleviates the “crowding problem”.

2.4. RF

An RF is an ensemble classifier consisting of a collection of tree structured classifiers $\left\{C\right(X,{\theta }_k)$, $k=$ 1, ...$\}$, where ${\theta }_k$ is defined as the independent identically distributed random vectors, and each decision tree casts a unit vote for the most popular class at input X.

A general RF framework is shown in Fig. 4 and is described as follows:

1) By employing bootstrap sampling, $k$ samples are selected from the training set, and the sample size of each selected sample is the same as those of the training sets.

2) $k$ decision tree models are built for $k$ samples, and $k$ classification results are obtained by these decision tree models.

3) According to $k$ classification results, the final classification result is decided by voting on each record. The RF increases the differences among the classification models by building different training sets; therefore, the extrapolation forecasting ability of the ensemble classification model is enhanced. After $k$ times training, a classification model series {$h_{1\left(X\right)}$, $h_{2\left(X\right)}$, …, $h_{k\left(X\right)}$} is obtained; the series is utilized to structure a multi-classification model system. The final classification result of the system is simple majority voting, and the final classification decision is given as Eq. (1):

where $H\left(x\right)$ is the ensemble classification model, hi is a single decision tree classification model, $Y$ is the objective output, and $I$ is an indicative function. Eq. 14 describes how the final classification is decided by majority voting.

Fig. 4Framework of the RF

3.1. Fault diagnosis for a hydraulic pump

The testing equipment for the axial piston hydraulic pump is shown in Fig. 5. In the experiment, the rotation speed is set at 5280 r/min, and the corresponding spindle frequency is 88 Hz. An accelerograph is installed at the end face of the pump. The sampling frequency is 1 kHz. The collected data contains three fault modes: normal, piston shoes and swashplate wear, and valve plate wear. The original signals of the three fault modes are shown in Fig. 6. Each data sample consists of 1024 points for analysis. Sixty samples are collected for each fault pattern. Thirty samples of each fault pattern are set as the training data, and the other 30 samples are set as the testing data.

Fig. 5A hydraulic pump system

Fig. 6Original signal of the hydraulic pump system

a)

b)

c)

The LCD method is applied to decompose the signals of the normal mode, piston shoes and swashplate wear mode, and valve plate wear mode. To acquire the fault feature vectors, the SampEn is utilized to quantify the ICSs. In this paper, the parameters are set to $m=$ 2 and $r=$ 0.2 times the standard deviation of the original signal.

Table 1 shows an example of the feature vectors of each fault mode represented by SampEn of ISCs. As shown in Table 1, because of the high dimension of ISC, the feature vector is usually composed of more than 8 SampEn values, which indicates the dimension of the feature vector is over 8 and makes the recognition of the fault mode difficult.

In this paper, t-SNE is used to reduce the dimension of the feature vector. Table 2 shows the results of applying t-SNE to the feature vectors. After using t-SNE, the feature vectors are reconstructed and reduced automatically to 3 dimensions. Through comparison of Table 1 and Table 2, it is obvious that the reconstructed feature vectors exhibit a strong ability of separability after dimension reduction by t-SNE, as shown in Fig. 7. Fig. 7 shows the results of dimension reduction using t-SNE. It provides a desirable input for the classification.

Fig. 7The feature vector extracted: a) using t-SNE, b) without t-SNE

a)

b)

Subsequently, the RF is applied to recognize the fault mode of the hydraulic pump. The results are shown in Table 3 and Fig. 8.

From the results diagnosed by the RF, we can see that the proposed method achieves the total accuracy of 98.89 %, in which the accuracy of the normal mode, the valve plate wear mode, and the piston shoes and swashplate wear mode is 96.67 %, 100.00%, and 100.00 %, respectively. However, the method without t-SNE only achieves the accuracy of 92.22 %, and the accuracy of the normal mode, the valve plate wear mode, and the piston shoes and swashplate wear mode is 90.00 %, 93.99 % and 93.33 %, respectively. From the results, it is obvious that after dimension reduction by t-SNE, the accuracy improved greatly. This fully demonstrates the effectiveness of the proposed method.

Fig. 8Comparison results of the diagnostic methods

3.2. Fault diagnosis for a bearing

The data of the bearing are from Xi’an Jiaotong University. In the experiment, the sampling frequency is 20 kHz. The data collected contains the following fault modes: normal, bearing inner race wear, bearing outer race wear and bearing rollers wear. The original signals are shown in Fig. 9. A data set of 1024 points for each team is selected for analysis. Then, 100 samples were collected for each fault pattern. Forty samples of each fault pattern are set as the training data, and the other 60 samples are set as the testing data.

Fig. 9Original signals of the bearing

a)

b)

c)

d)

Table 4 presents the obtained SampEn values of these ISCs.

Subsequently, the original features of the datasets are reduced automatically to 3 dimensions by t-SNE. For comparison, the dimension of the feature vector is also reduced by PCA and LPP, as shown in Fig. 10. From Fig. 10, we can conclude that the feature vector reduced by t-SNE exhibits a stronger ability for separability than those reduced by PCA and LPP and thus provides a desirable input for classification.

Subsequently, the RF was applied to recognize the states of the hydraulic pump based on the feature vectors obtained by LCD-SampEn and t-SNE. The result is shown in Table 5 and Fig. 11. A comparison of the results obtained by the PCA and LPP methods is also shown in Table 5.

Fig. 10The feature vector extracted using: a) PCA, b) LPP, c) t-SNE

a)

b)

c)

From the results diagnosed by the RF, we can see that the PCA method shows the worst performance of the three methods; the accuracy rates of the normal mode, the inner race wear mode, the outer race wear mode, and the rollers wear mode are 86.67 %, 90.00 %, 90.00 % and 91.67 %, respectively, and the total accuracy is 89.59 %. The poor results are due to the non-linear nature of the feature, which is not ideal for use with the PCA method, which is a typical linear dimension reduction method. For the LPP method, the accuracy rates of the normal mode, the inner race wear mode, the outer race wear mode, and the rollers wear mode are 91.67 %, 90.00 %, 93.33 % and 91. 67 % respectively, and the average accuracy is 91.67 %. The t-SNE method shows the best performance of the three methods; the accuracy rates of the normal mode, the inner race wear mode, the outer race wear mode, and the rollers wear mode are 98.33 %, 100.00 %, 98.33 % and 98.33 %, respectively, and the total accuracy is 98.75 %. The result verifies the superiority of the t-SNE method in feature dimension reduction and fault diagnosis.

Fig. 11Comparison results of the diagnostic methods