Diagnosis of rolling element bearing fault arising in gearbox based on sparse morphological component analysis

1. Introduction

It is meaningful to extract the fault feature of rolling element bearing effectively because the rolling bearing is one of the most commonly used and easily damaged mechanical components. In recent years, various signal processing techniques have been developed to extract the fault features of machinery. Traditional statistical features on time series are usually used to monitor machinery condition [1, 2]. Frequency domain method can identify frequency structure of signal, and it’s widely used in gearbox signal analysis [3, 4]. After analyzing the shortcomings of current feature extraction and fault diagnosis technologies of rolling bearing, the wavelet packet decomposition and empirical mode decomposition were combined in paper [5] to extract fault feature frequency, and the extracted feature was used as the input of neural network. At last the rotating machinery’s early fault was diagnosed successfully. A new bearing fault diagnosis method named cyclic spike detection method was proposed to extract the features of bearing’ weak fault [6]. Two multi-component signal mixtures, including a simulated signal and a real vibration signal collected from an industrial machine were used to validate the effectiveness of the proposed method. The results demonstrated that the proposed method can extract the bearing weak fault features from other strong masking vibration signals and noise. The minimum variance cepstrum (MVC) was used to detect the early fault in automotive ball bearing and the experiment results demonstrated that MVC could determine bearing fault periods more clearly than other methods under the given condition [7]. The discrete wavelet transform (DWT) and artificial neural network (ANN) were combined to diagnose gear box incipient fault in paper [8]. The extracted features from DWT were used as the inputs of a neural network for classification purposes. The results showed that the developed method can diagnose different conditions of the gear box reliably. Q. Hai et al. compared the performance of wavelet decomposition-based de-noising method and wavelet filter-based de-noising method based on signals collected from mechanical defects [9]. The comparison results revealed that wavelet filter-based de-noising method is much more suitable and more reliable to detect a weak signature of mechanical impulse-like defect signals, whereas the wavelet decomposition de-noising method can achieve satisfactory results on smooth signal detection. Aiming at the difficulty of early weak impulse fault feature extraction of rolling bearings, an integrated algorithm based on resonance demodulation method and entropy threshold de-noising of wavelet packet coefficients was proposed to determine the resonance modulation bands resulted by early weak impulse faults signals accurately [10]. The validity and practicability of the proposed method in feature extraction of rolling bearing’ early weak impulse fault were proved through experiment.

MCA based on sparse decomposition uses basic pursuit to separate the signal components by utilizing multiple dictionaries [11]. The assumption about MCA lies that each dictionary is efficient in sparsely representing one component but being highly inefficient in representing others. MCA has many successful applications on image processing [12, 13] and source separation [14, 15]. Vibration signal usually contains complex components from different sources, and these components are always morphologically different. By separating components from vibration signal using MCA method, fault related signal can be distinguished from the rest parts.

The paper is organized as follows. Section 2 is dedicated to sparse representation and MCA. Section 3 is simulation verifying the effectiveness of the proposed method. In Section 4 the experimental analyzed results of rolling element bearing fault arising in gearbox are presented. In Section 5, the analyzed experimental results using independent component analysis (ICA) as comparison are presented and the advantages of the proposed method are verified. Conclusions obtained from the above results are given in Section 6.

2. Morphological component analysis based on sparse representation

2.3. The construction of over-complete dictionaries for sparse decomposition of gear and bearing fault signals

According to the dynamic model of fixed-axis gear system, its vibration equation can be expressed by [16]:

3

$\left\{\begin{array}{l}{x}_g\left(t\right)=\sum _{m=1}^M{X}_m\left[1+a_m\left(t\right)\right]\mathrm{c}\mathrm{o}\mathrm{s}\left[m2\pi f_{z}t+b_m\left(t\right)+{\phi }_m\right]+n\left(t\right),\\ a_m\left(t\right)=\sum _{n=1}^N{A}_{m,n}\mathrm{s}\mathrm{i}\mathrm{n}\left[n2\pi f_{r}t+{\alpha }_{m,n}\right],\\ b_m\left(t\right)=\sum _{n=1}^N{B}_{m,n}\mathrm{s}\mathrm{i}\mathrm{n}\left[n2\pi f_{r}t+{\beta }_{m,n}\right],\end{array}\right.$

where $x_g$ is the gear vibration, $f_z$ is the meshing frequency, $X_m$ is the amplitude of vibration, ${\phi }_m$ is the phase. $a_m\left(t\right)$ and $b_m\left(t\right)$ is the modulation function, $n\left(t\right)$ is the random noise. $A_{m,n}$ and $B_{m,n}$ are the amplitude and phase modulation factors. $f_r$ is the rotating frequency. $m$ is the order. ${\alpha }_{m,n}$ and ${\beta }_{m,n}$ are the initial phase angles.

Set the parameters in the gear model, and the time-domain and frequency-domain waveforms are shown in Fig. 1(a) and (b). Based on Fig. 1, the main components of gear vibration signal are cosines. So, the cosine packet over-complete dictionary is used for analysis of gear vibration signal.

Fig. 1Simulated signal of fault gear

a) Time-domain waveform

b) Frequency-domain waveform

Using the rolling bearing fault mode whose mathematical model [17, 18] is expressed in equation (4). ${\tau }_i$ is the tiny fluctuation around mean period $T$:

4

$\left\{\begin{array}{l}{x}_b\left(t\right)=s\left(t\right)+n\left(t\right)=\sum _i{A}_{i}h\left(t-iT-{\tau }_i\right)+n\left(t\right),\\ A_i=A_0\mathrm{c}\mathrm{o}\mathrm{s}\left(2\pi f_{r}t+{\varphi }_A\right)+C_A,\\ h\left(t\right)=e^{-Bt}\mathrm{c}\mathrm{o}\mathrm{s}\left(2\pi f_{n}t+{\varphi }_{\omega }\right).\end{array}\right.$

Set the parameters in the bearing fault model, and the time-domain and frequency-domain waveforms are shown in Fig. 2(a) and (b). Based on Fig. 2, the main components of bearing fault vibration signal are impulses. The symmlet8 wavelet packet is taking on impulse characteristic, so the symmlet8 wavelet packet dictionary is a good candidate for bearing fault signal impulse extraction.

Fig. 2Simulated signal of rolling element bearing

a) Time-domain waveform

b) Frequency-domain waveform

2.4. MCA using multiple dictionaries

Morphological component analysis can decompose a signal into its building parts. Assume that the signal $s$ is a linear combination of $K$ components ${\varphi }_k$, $s=\sum _{k=1}^K{\varphi }_k$, where ${\varphi }_k$ represents each type of signals to be decomposed. The MCA model also assumes that a possibly over-complete dictionary ${\mathrm{\Phi }}_k$ exists for each $k$, in which each component admits a sparse representation while its decomposition over other ${\mathrm{\Phi }}_{k\text{'}\ne k}$ is non-sparse. Ideally, the MCA model can be described as follows:

5

$\underset{\{{\alpha }_1,\dots ,{\alpha }_K\}}{\mathrm{m}\mathrm{i}\mathrm{n}}\sum _{k=1}^K{‖{\alpha }_k‖}_0\mathrm{}\mathrm{}\text{subject to}\mathrm{}s=\sum _{k=1}^K{\mathrm{\Phi }}_k{\alpha }_k.$

This optimization task is likely to lead to a successful separation of the signal. This expectation relies on the assumption made earlier about the over-complete dictionary ${\mathrm{\Phi }}_k$ being efficient in sparsely representing ${\varphi }_k$ and being highly inefficient in representing other components.

The minimization of $l_0$ norm in Eq. (5) is still a non-convex optimization problem and difficult to solve. So, the $l_1$ norm is used in the optimization problem. In addition, due to additive noise or model dis-match, signal will generally not decompose cleanly into $K$ components. We propose the alternative optimization problem by relaxing the constraint as follows:

6

$\underset{\{{\alpha }_1,\dots ,{\alpha }_K\}}{\mathrm{m}\mathrm{i}\mathrm{n}}\sum _{k=1}^K{‖{\alpha }_k‖}_1+\lambda {‖s-\sum _{k=1}^K{\mathrm{\Phi }}_k{\alpha }_k‖}_2^2.$

The signal can be decomposed into signal-plus-residual with sparse representation:

7

$s=f^{\left(\lambda \right)}+r^{\left(\lambda \right)},$

where $f^{\left(\lambda \right)}=\sum _{k=1}^K{\mathrm{\Phi }}_k{{\alpha }_k}^{\left(\lambda \right)}$, and $r^{\left(\lambda \right)}=s-\sum _{k=1}^K{\mathrm{\Phi }}_k{\alpha }_k$. The residual is controlled by the parameter $\lambda$. For an exact representation of the data, $\lambda$ must be set to zero. But when noise is contained in data, $\lambda$ should be set to few times the noise standard deviation $\sigma$. Assuming the dictionary is normalized, we usually set $\lambda =\sigma \sqrt{2\mathrm{l}\mathrm{o}\mathrm{g}\left(p\right)}$, where $p$ is the cardinality of the dictionary.

Considering that the dictionary is over-complete, storing and manipulating coefficients ${\alpha }_k$ requires large memory. We can obtain a simplified optimization problem:

8

$\underset{\{s_1,\dots ,s_K\}}{\mathrm{m}\mathrm{i}\mathrm{n}}\sum _{k=1}^K{‖{\mathrm{\Phi }}_k^{+}{s}_k‖}_1+\lambda {‖s-\sum _{k=1}^K{s}_k‖}_2^2,$

where $s_k={\mathrm{\Phi }}_k{\alpha }_k$; ${\mathrm{\Phi }}_k^{+}$ is Moore-Penrose pseudo-inverse matrix. The optimization problem above can be used to solve each component directly, which helps to reduce the computational complexity and save the computer memory.

The optimization problem in Eq. (8) is still complex to solve. By using Block Coordinate Relaxation (BCR) method [18], the above optimization problem can be replaced by a sequence of easier problems. Assuming that $\{{\alpha }_k{\}}_{k=2}^K$ are all known, we can get a simplified task as:

9

$\underset{\{{\alpha }_1,\dots ,{\alpha }_K\}}{\mathrm{m}\mathrm{i}\mathrm{n}}{‖{\alpha }_1‖}_1+\lambda {‖\left(s-\sum _{k=2}^K{\mathrm{\Phi }}_k{\alpha }_k\right)-{\mathrm{\Phi }}_1{\alpha }_1‖}_2^2.$

Here, an iterative thresholding algorithm is proposed to update coefficients, using a threshold which decreases linearly towards zero with iterations. As displayed in Fig. 3, the MCA algorithm based on BCR is as follows:

(1) Initialize the number of iterations $L_{\mathrm{m}\mathrm{a}\mathrm{x}}$ and threshold $\delta =\lambda ·L_{\mathrm{m}\mathrm{a}\mathrm{x}}$.

(2) Perform $J$ times:

For $k=\mathrm{1,2},\dots ,K$:

– Calculate the residual $r_k=s-{\sum }_{k\text{'}\ne k}{s}_{k\text{'}}$;

– Calculate the coefficients ${\alpha }_k={\mathrm{\Phi }}_k^{+}{r}_k$;

– Threshold the coefficients ${\alpha }_k$with the threshold $\delta$and obtain ${\widehat{\alpha }}_k$;

– Reconstruct component $s_k$ by $s_k={\mathrm{\Phi }}_k{\widehat{\alpha }}_k$.

(3) Update the threshold exponentially.

(4) If $\delta >\lambda$, return to Step (2). Else, finish.

In Step (2), a hard threshold towards the end of the iterative process may lead to better results compared with soft threshold, which is given by:

10

${\widehat{\alpha }}_k=\left\{\begin{array}{ll}{\alpha }_k,& \left|{\alpha }_k\right|\ge \delta ,\\ 0& \left|{\alpha }_k\right|<\delta .\end{array}\right.$

Fig. 3The flowchart of MCA algorithm based on BCR

3. Simulation

The study focus of the paper is the performance of MCA in separating the collected signal from gearbox which containing both gear vibration and fault rolling element bearing vibration. The simulation vibration signal can be represented by the time-domain addition of Eq. (3) and Eq. (4) as following:

Set the sampling frequency $f_s$ as 10 kHz, and the parameters of rolling element bearing inner fault are set as follows: the rotating frequency $f_r=$50 Hz, the natural frequency $f_n=$ 2 kHz and the inner race frequency $f_o=$96 Hz. The parameters of gear vibration are set as follows: the moderating frequencies of the amplitude moderated and frequency moderated are $f_r$, 2$f_r$, 3$f_r$, 4$f_r$ and 5$f_r$, and the corresponding moderating coefficients are 1, 0.8, 0.7, 0.6, 0.5.

The time-domain waveforms of the gear vibration and rolling element bearing inner fault vibration are shown in Fig. 4(a) and (b). Fig. 4(c) is the random white noise, and the combined signal of (a), (b) and (c) is shown in Fig. 4(d).

Envelope analysis based on Hilbert transformation is used for demodulation, which is usually used to extract fault characteristic frequencies in vibration signal.

Hilbert transformation $h\left(t\right)=H\left\{x\left(t\right)\right\}$ of a real signal $x\left(t\right)$ is defined as:

Fig. 4Simulated signal of fault gearbox

a) Time-domain waveform of gear

b) Time-domain waveform of rolling element bearing inner race fault

c) Random white noise

d) Combined signal

The analytical signal $z\left(t\right)$ is defined as the combination of $x\left(t\right)$ and its Hilbert signal as:

The envelope $E\left(t\right)$of the given signal $x\left(t\right)$ is defined as:

By performing FFT on the envelope $E\left(t\right)$, we can obtain the envelope spectrum of the given signal. Fault characteristic frequencies in this envelope spectrum can indicate machinery fault features.

Apply envelope demodulation spectrum to the signal shown in Fig. 4(d) and the analyzed result is given in Fig. 5 from which the inner race fault characteristic frequency could not be extracted.

Fig. 5The envelope demodulation spectrum of the signal shown in Fig. 4(d)

Fourier dictionary is effective in analyzing harmonic components, and the symmlet8 wavelet packet atoms can match the bearing fault model well as the above stated. So, these two types of dictionaries are exploited to extract gear vibration and bearing fault vibration.

Two components are extracted from the raw signal after using MCA method. The first component in Fig. 6(a) which is reconstructed with the Fourier dictionary. It can be defined as the gear vibration signal. Its PSD estimation gives a detailed frequency structure in Fig. 6(b). We can clearly see the gear meshing frequency and its sidebands. In Fig. 6(c), the envelope demodulation spectrum of the extracted signal shown in Fig. 6(a) is given, and the main components of spectrum structure are rotating frequency with its harmonics.

Fig. 6Gear signal separated by MCA

a) The separated gear signal using MCA

b) The fft of the signal shown in Fig. 6(a)

c) The envelope demodulation spectrum of the signal shown in Fig. 6(a)

The second separated component in Fig. 7(a) which is reconstructed with the Symmlet8 wavelet packet dictionary. It identifies transient features from rolling element bearing inner race fault. Its fft and envelope analysis are presented in Fig. 7(b) and (c). The inner race fault feature frequency and its harmonics are extracted successfully based on Fig. 7(c). Besides, the moderating frequency (rotating frequency) is also extracted. Much better separated results will be obtained if updating the threshold exponentially as shown the step 3 in Fig. 3 [19], and the Fig. 8 shows the threshold value decreases exponentially. Fig. 9 shows that the residual error decreases quickly with the number of iteration.

Fig. 7Bearing inner race fault signal separated by MCA

a) The separated rolling element bearing inner race fault signal using MCA

b) The fft of the signal shown in Fig. 7(a)

c) The envelope demodulation spectrum of the signal shown in Fig. 7(a)

Fig. 8Threshold updating with an exponential method

Fig. 9The residual error for reconstructed signal using MCA

The simulation example shows that different components can be extracted successfully from vibration signal collected for gearbox by selecting proper dictionaries, which offers a potential way to detect machinery fault arising in gearbox.

4. Experiment

The experimental gearbox test rig is shown in Fig. 10. During the test, the motor speed is adjusted by the inverter control. The magnetic power brake is used to adjust the load. The gearbox is a single stage transmission, and the two shafts of the gear box are supported by four rolling bearings of the same type of NSK#1302.There are four acceleration sensors and one data acquisition card USB9234 in the data acquisition system. The programming of data acquisition software is realized by National Instruments LabView. The transmission structure and the arrangement of measuring points are shown in Fig. 11.

The inner race of the 2nd measuring point bearing is processed pitting failure using EDM (Electron Discharge Machining) in the experiment, and the other three measuring points of the bearings are normal. The pitting failure is located in the inner race of the test bearing and the size of the failure is 1 mm×1 mm. The sampling frequency is 44.1 kHz, and the rotation frequency $f_r$ of the input shaft is 100/6 Hz. As can be seen from Fig. 1, the transmission ratio is $i=Z_2/Z_1=$39/28, and the gear modulus is 2. The rated power of variable frequency of motor is 0.55 kw. The geometric parameters of the fault bearing are shown in Table 2, and the relationship of the theoretical characteristic frequencies of the fault bearing and the rotation frequency $f_r$ are shown in Table 3.

Fig. 10The test rig

Fig. 11The structure of the transmissions and the layout of measuring points

The time-domain waveform and the envelope spectrum of the fault signal are shown in Fig. 12(a) Fig. 12(b). We could not identify the bearing fault basing on Fig. 12(a) because the impulse characteristic is buried in the strong background noise. In the Fig. 12(b), the spectrum lines of the spectrum are very chaotic, which could not reflect the fault characteristic information. By the above analysis results, although the acceleration sensor is so close to the fault bearing, using the traditional signal processing methods could not effectively extract the bearing fault signal characteristics. The main reasons can be summarized as follows: the vibration signal of the gear box is complex which contains the gear meshing signal components, the vibration of shaft, the vibration of fault bearing and the vibration of the motor, etc. Therefore, it is very necessary and significant to further separate the signal of the gear box in order to better extract the fault characteristic signals.

The analysis results of the experimental signal shown in Fig. 12(a) based on MCA are given. Same as the discussions in Section 2 and 3, the Cosine packet over-complete dictionary and Symmlet8 wavelet packet over-complete dictionary are used to match the gear vibration signal and bearing fault vibration signal separately. Besides, the updating of the threshold exponentially is also same as the Fig. 8. Two components are obtained after MCA analysis and their corresponding time-domain waveforms are shown in Fig. 12(c) and (e). The second component shown in Fig. 12(e) can be identified as the bearing fault vibration signal because it takes on impulse characteristic, and the characteristics of the second component shown in Fig. 12(c) is much close to the gear vibration signal.

Fig. 12MCA analysis of the gearbox signal with an inner race defect bearing

a) The time-domain waveform of the gearbox signal with an inner race fault bearing

b) The envelope demodulation spectrum of the signal shown in Fig. 12(a)

c) Component 1 separated using MCA

d) The envelope demodulation spectrum of the signal shown in Fig. 12(c)

e) Component 2 separated using MCA

f) The envelope demodulation spectrum of the signal shown in Fig. 12(e)

The envelope demodulation spectrum analysis results of the two components are also given in Fig. 12(d) and (f) respectively in order to reflect the separated results much clearer. In Fig. 12(f) the values of spectral lines are mainly distributing on 105 Hz and its harmonics. It can judge the 105 Hz is bearing inner race fault characteristic frequency through comparison the fault characteristic frequencies shown in Table 2 (There exists 1.1 Hz error due to the manufacture and installation error). The effectiveness of MCA in separating gear vibration and bearing fault vibration signal is verified again through the above experimental signal analysis results.

Same as the above experiment and analysis process, a bearing with an outer race fault is installed near the 2nd measuring point, the other three bearings are normal, and the corresponding analyzed results are given in Fig. 13. Though the envelope demodulation analyzed shown in Fig. 13(b) of the original fault signal shown in Fig. 13(a) can extract the outer race characteristic frequency roughly, the analyzed result shown in Fig. 13(f) is much better after MCA handling.

Fig. 13MCA analysis of the gearbox signal with an outer race defect bearing

a) The time-domain waveform of the gearbox signal with an outer race fault bearing

b) The envelope demodulation spectrum of the signal shown in Fig. 13(a)

c) Component 1 separated using MCA

d) The envelope demodulation spectrum of the signal shown in Fig. 13(c)

e) Component 2 separated using MCA

f) The envelope demodulation spectrum of the signal shown in Fig. 13(e)

5. Comparison

Independent component analysis is a well-known method for the blind source separation while assuming that the sources are statistically independent and non-Gaussian [20]. But ICA is unsuitable for underdetermined cases where sensors are less than sources. Therefore, it cannot separate multiple sources from only single-channel data. Here, the analysis results of the gearbox signal with an inner race fault bearing using ICA are given in Fig. 14 as comparison 1. Neither of them can detect bearing fault feature in gearbox signal.

The analysis results of the gearbox signal with an inner race fault bearing using MP are given in Fig. 15 as comparison 2. Neither of them can detect bearing fault feature in gearbox signal same as the analysis results using ICA method.

Fig. 14Analysis result using ICA

a) Component I

b) Envelope analysis of Component I

c) Component II

d) Envelope analysis of Component II

Fig. 15Analysis result using MP

a) Component I

b) Envelope analysis of Component I

c) Component II

d) Envelope analysis of Component II

6. Conclusions

The principle and algorithm of MCA method basing on sparse representation theory is discussed in the paper, then the characteristics of proposed method is described by the simulation, and the separation different morphological constituents ability of the proposed method is verified through simulation and experiment. Besides, compared with other transformations such as FFT or wavelet, the MCA can represent a signal with much lesser atoms by using proper over-complete dictionary, which can capture the underlying signal structure and save storage.

With morphological diversity in dictionaries, MCA can separate a given signal into its additive components. It solves single channel blind source separation problem that ICA cannot deal with, and the comparison result of the proposed method with ICA is given at last. The advantage of the proposed method is verified furthermore.