Application of EMD-AR and MTS for hydraulic pump fault diagnosis

2.1. EMD method

Huang et al. [5] developed EMD in 1998. The EMD method is developed based on the simple assumption that any signal consists of different simple intrinsic modes of oscillations. Each linear or non-linear mode will have the same number of extrema and zero-crossings. Only one extremum exists between successive zero-crossings. Each mode should be independent of the others. In this way, each signal could be decomposed into a number of IMFs, each of which must satisfy the following definitions [9]: (1) In the whole data set, the number of extrema and the number of zero-crossings must either be equal or differ by at most one; (2) At any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero.

An IMF represents a simple oscillatory mode compared with the simple harmonic function. With the definition, any signal $x\left(t\right)$ can be decomposed as follows:

Step 1: Identify all the local extrema, and then connect all the local maxima by a cubic spline line to give the upper envelope.

Step 2: Repeat this procedure for the local minima to produce the lower envelope. Between them, the upper and lower envelopes should cover all the data.

Step 3: The mean of the upper and lower envelope values is designated as $m_1\left(t\right)$; the difference between the signal $x\left(t\right)$ and $m_1\left(t\right)$ is the first component $h_1\left(t\right)$:

Ideally, if $h_1\left(t\right)$ is an IMF, then $h_1\left(t\right)$ is the first component of $x\left(t\right)$.

An equivalent set of two first-order non-autonomous equations is as follows:

Step 4: If $h_1\left(t\right)$ is not an IMF, $h_1\left(t\right)$ is treated as the original signal and steps 1, 2 and 3 are repeated; then:

After repeated sifting, i.e., up to $k$ times, $h_{1k}\left(t\right)$ becomes an IMF, that is:

Then, it is designated as:

The first IMF component from the original data $c_1\left(t\right)$ should contain the finest scale or the shortest period component of the signal.

Step 5: Separating $c_1\left(t\right)$ from $x\left(t\right)$, we get:

$r_1\left(t\right)$ is treated as the original data, and the above processes are repeated; therefore, the second IMF component $c_2\left(t\right)$ of $x\left(t\right)$ can be obtained. Repeating the process described above $n$ times, $n$-IMFs of signal $x\left(t\right)$ can be obtained. Then:

The decomposition process can be stopped when $r_n\left(t\right)$ becomes a monotonic function from which no more IMFs can be extracted.

Using this procedure, any signal can be decomposed. We finally obtain:

Thus, the signal is decomposed into $n$-empirical modes and a residue $r_n\left(t\right)$, which is the mean trend of $x\left(t\right)$. Each IMF $c_1\left(t\right),c_2\left(t\right),\dots ,\mathrm{}{c}_n\left(t\right)$ contains lower-frequency oscillations than the prior-extracted one, while $r_n\left(t\right)$ represents the central tendency of signal $x\left(t\right)$.

EMD has been shown to be a fast, effective, self-adaptive method for nonlinear and non-stationary time series analysis, but there is still a drawback worth noting: the end effect, whereby distortion appears at the end of the signal in the decomposition process. In this paper, we employ the mirror periodic extending method (MPM) to solve this problem.

2.2. MTS

The MTS starts by collecting data on normal observations. The MD is calculated using certain characteristics to determine whether the MD has the ability to differentiate a normal group from an abnormal group. If the MD cannot identify the normal group using those particular characteristics, then a new combination of characteristics are needed to be examined. When the correct set of characteristics is determined, the Taguchi methods are employed to evaluate the effect of each characteristic. If possible, dimensionality is reduced by eliminating those characteristics that do not add value to the analysis. The MTS consists of the following three stages [8, 10].

2.3. Diagnosis approach for hydraulic pump

First, the fault vibration signal of the hydraulic pump $X_t$, is decomposed by EMD into $n$ IMFs, $im{f}_1,im{f}_2,...,im{f}_n$. Each component represents a different characteristic information. After the EMD method is applied to $X_t$, the IMFs can completely combine the characteristics of $X_t$Therefore, the characteristics of $X_t$can be obtained by extracting the characteristics of $im{f}_1,im{f}_2,...,im{f}_n$ [9].

Second, the following AR model, $AR\left(m\right)$, is established for each IMF component as:

where ${\phi }_{ik}\mathrm{}(k=\mathrm{1,2},...,m)$ represent the parameters of the model, AR indicates an AR model of order $m$, $c_i\left(t\right)$ denotes the remnant of the model and is a white-noise sequence whose mean value is zero and variance is ${\sigma }_i^2$. The parameters ${\phi }_{ik}\left(k=\mathrm{1,2},...,m\right)$ can reflect the inherent characteristics of a hydraulic pump vibrating system. The variance of the remnant ${\sigma }_i^2$ is closely related with the output characteristics of the system. ${\phi }_{ik}(k=\mathrm{1,2},...,m)$ and ${\sigma }_i^2$ can be chosen as the initial feature vectors $A_i=[{\phi }_{i1},{\phi }_{i2},...,{\phi }_{im},{\sigma }_i^2]$.

The singular values of the feature vector matrix are calculated to extract the feature vectors. Based on matrix theory, the singular value of the matrix is the inherent characteristic of the matrix and has good stability, that is, when the matrix elements encounter small changes, the matrix singular value only slightly changes. Thus, the singular values of the feature vector matrix obtained by AR can be extracted as mechanical component characteristics. However, the SVD also has its disadvantages, that is, for a space reconstruction in the time series, embedding dimension and delay decisions have no specific theoretical basis. The initial feature vector matrix obtained by AR could avoid this problem. Accordingly, the initial feature vector $A_i$, from $n$ IMF components constitutes the initial feature vector matrix $A$, where $A={(A_1,A_2,...,A_N)}^T$. $A$ is decomposed by SVD to extract the singular value ${\sigma }_1,{\sigma }_2,...,{\sigma }_N$, thereby extracting the fault feature vectors. The criterion of condition identification is the MD. The proposed diagnosis method is illustrated in Fig. 1.

The fault diagnosis method for the hydraulic pump is as follows:

(1) Acquire $n$ vibration signal samples at a certain sampling frequency $f_s$, under each of the following conditions: hydraulic pump is normal, with valve plate wear, and with slipper loosing. The $3n$ signals are taken as the samples.

(2) Each signal is decomposed by EMD, and finite number of IMF components can be obtained.

(3) Normalize each IMF component to obtain a new component to eliminate the effect of the signal amplitude on the variance of the remnant ${\sigma }_i^2$:

Fig. 1Flow chart of the proposed method

(4) Construct the AR model for the normalized component. Determine the order $m$, of the model using the final prediction error (FPE) criterion and estimate the AR parameters ${\phi }_{ik}\left(k=\mathrm{1,2},\dots ,m\right)$ and the remnant’s variance ${\sigma }_i^2$ by the minimum squares method, where ${\phi }_{ik}$ denotes the $k^{th}$ AR parameters of the $i^{th}$ IMF component. ${\phi }_{ik}\left(k=\mathrm{1,2},\dots ,m\right)$ and ${\sigma }_i^2$can be determined. These values can be combined to construct the initial feature vector of the $i^{th}$ IMF component as follows:

where $j=$1, 2, 3 denotes the normal condition, condition with valve plate wear, and condition with slipper loosing.

(5) Each sample in all conditions is decomposed by SVD to obtain $m+1$ singular values. Formally, the SVD of an $m\times n$ real or complex matrix $\mathbf{M}$ is a factorization of the form:

where $U$ is an $m\times m$ real or complex unitary matrix, $\sum$ is an $m\times n$ rectangular diagonal matrix with nonnegative real numbers on the diagonal, and ${\mathbf{V}}^{\mathrm{*}}$ (the onjugate transpose of $\mathbf{V}$) is an $n\times n$ real or complex unitary matrix. The diagonal entries ${\sum }_{i,i}$ of $\sum$ are known as the singular values of $\mathbf{M}$.

(6) Determine the MD between the testing data and the benchmark space composed of the normal training data as follows:

where ${\mathbf{R}}^{-1}$ is the inverse of the correlation matrix $\mathbf{R}$, $M{D}_j$ is the MD for the $j^{th}$ observation, and $k$ is the number of characteristics.

(7) Establish the orthogonal table, optimize the benchmark space, reduce the redundant features, and extract the principal characteristics by Taguchi methods. The MD values for all of the datasets are re-calculated by using the optimized algorithm.

(8) Identify the fault condition of the hydraulic pump.

3.1. Experimental setup

In this study, a test plunger pump rig, as shown in Fig. 2, was tested and analyzed to verify the presented fault diagnosis method. In the experiment, two commonly occurring faults in the plunger pump were set, namely, slipper loosing and valve plate wear. Under three conditions, including the two faulty conditions and the normal state, the vibration signal was acquired from the end face of the plunger pump with a stabilized motor speed of 528 r/min, and sampling rate of 1000 Hz. Twelve samples were acquired under the normal condition and four samples were obtained for each of the two fault states. Among these samples, the first eight normal samples were used to construct the benchmark space, and the remaining samples were used for testing.

Fig. 2Experimental test plunger pump rig

3.2. Feature extraction

The feature vectors were determined by the proposed method. The vectors were calculated using the first eight IMF components. Fig. 3 shows the acceleration vibration signal of the hydraulic pump with a normal signal. It is decomposed into eight IMFs by EMD. As shown in Fig. 4, each IMF component has a distinct time characteristic scale.

The system condition is mainly decided by the first several AR parameters and the remnant variance. The order of the model $m$, was determined by FPE criterion. This order is different from the number of IMF components, wherein the maximum and minimum components are 26 and 9, respectively. In this case, the first seven AR parameters ${\phi }_{ik}\left(k=1,2,\dots ,7\right)$ and ${\sigma }_i^2$ were chosen. These values constitute the eight-dimension vector. The extracted feature vectors are listed in Table 1 (only the feature vectors of the first four IMF components from one sample in each condition are listed in Table 1 because of space limitations). ${\mathbf{A}}_{ji}\mathit{}\left(j=1,2,3;i=1,2,3,\dots ,8\right)$ denotes the feature vector of the $i^{th}$ IMF component under the $j^{th}$ condition, where $j$ denotes the normal condition, condition with valve plate wear, and condition with slipper loosing, and $i$ denotes the first eight IMF components. The feature vector matrixes were decomposed to obtain the singular value by SVD. The extracted singular values are listed in Tables 2 and 3. The extracted singular values for the model training (Table 2) are used to construct the benchmark space in calculating the MD, and the others (Table 3) are used for testing.

Fig. 3Acceleration vibration signal of the hydraulic pump with a normal signal

Fig. 4Decomposed results of the hydraulic pump vibration signal shown in Fig. 3 by EMD

3.3. Pattern recognition