Abstract

The generated signals generally contain a large amount of background noise when the mechanical bearing fails, and the fault signals present nonlinear and non-Gaussian feature, which have heavy tail and belong to -stable distribution (); even the background noises are also -stable distribution process. Then it is difficult to obtain reliable conclusion by using the traditional bispectral analysis method under -stable distribution environment. Two improved bispectrum methods are proposed based on fractional lower-order covariation in this paper, including fractional low-order direct bispectrum (FLODB) method, fractional low-order indirect bispectrum (FLOIDB) method. In order to decrease the estimate variance and increase the bispectral flatness, the fractional lower-order autoregression (FLOAR) model bispectrum and fractional lower-order autoregressive moving average (FLOARMA) model bispectrum methods are presented, and their calculation steps are summarized. We compare the improved bispectrum methods with the conventional methods employing second-order statistics in Gaussian and distribution environments; the simulation results show that the improved bispectrum methods have performance advantages compared to the traditional methods. Finally, we use the improved methods to estimate the bispectrum of the normal and outer race fault signal; the result indicates that they are feasible and effective for fault diagnosis.

1. Introduction

Bispectral analysis based on high-order statistics is an effective tool to solve nonlinear phase coupling and non-Gaussian fault diagnosis [1, 2]. The traditional bispectrum methods include nonparametric bispectrum [3, 4], parametric AR bispectrum [5, 6], parametric ARMA bispectrum [7], and their improved bispectrum methods [8]. The bispectrum of the signal contains not only the amplitude information but also the phase information. Bispectrum can effectively suppress the influence of Gaussian background noise and extract the non-Gaussian features hidden in the signal, and the graphics are intuitive. The fault feature of different nonlinear coupling modes in the bispectrum can be applied to quickly identify the working state of the bearing. Therefore, bispectrum can better extract the signal features than the traditional power spectrum, which has been widely used in mechanical fault diagnosis [911]. In recent years, the methods of mechanical fault signal analysis have been developed. A new family of model-based impulsive wavelets and their sparse representation method are presented for rolling bearing fault diagnosis in [12]. An SVD principle analysis method based on the correlation coefficient is proposed for the bearing fault diagnosis in [13]. Subsequently, Qin et al. proposed a K-SVD algorithm with adaptive transient dictionary and transient feature extraction by the improved orthogonal matching pursuit [14]. Guo et al. applied the resonance demodulation and vibration separation to tooth root crack detection of planet and sun gears [15].

Recently, some extension methods based on the traditional bispectrum have also been studied and applied to the rotating machinery signal analysis, such as the deterministic bispectrum [16], modulation signal bispectrum [17, 18], principal component bispectrum [19, 20], and cyclic bispectrum [21]. Cheng et al. proposed a mechanical fault location and diagnosis method based on two bispectra and fuzzy clustering, which can effectively diagnose the fault state and location [22]. A new fault diagnosis method based on wavelet packet decomposition and modulation signal bispectrum analysis was proposed in [23]. The method firstly reconstructs the wavelet packet energy signal in time-frequency domain and then carries out modulation bispectrum analysis on the reconstructed signal, which can realize early fault diagnosis. Wang et al. proposed a bispectrum image texture features manifold method based on the support vector machines and genetic optimization algorithms for the rolling bearing vibration signal analysis [24]. A new bispectrum analysis method based on the optimal scale shape slice was proposed in [25], and the method can extract the fault signal from the sensitive modal component; hence, it can better extract fault features. A new variant modulation signal bispectrum method was introduced in [26], which was used to measure and analyze the fault current signals of the different mechanical motors, and the results show that this method is better than the traditional bispectrum method. The traditional and improved bispectrum methods have been applied in the rotating machinery signal analysis fields, but the methods still have some defects, and their performances degrade in impulsive environment and even fail. Therefore, it is of great significance to explore high efficiency and performance of bispectrum analysis methods.

In actual working conditions, the mechanical bearings work in poor environment, the generated signals generally contain a large amount of background noise when the fault occurs, and the fault signals have obvious nonlinear and non-Gaussian properties, which belong to -stable distribution process, and even the same with the noise in the signals [2730]. Hence, it is difficult to find a solid conclusion by using the traditional bispectrum analysis methods. Therefore, improved bispectrum methods which can be applicable to -stable distribution environment need to be explored. Recently, -stable distribution model was used for statistical modeling of the ocean environmental noise [31].

The adaptive cumulative distribution detector and blind estimation of frequency hopping parameters methods were proposed based on -stable distribution model in [32, 33]. Several improved frequency spectrum analysis methods have been introduced for -stable distribution environment in [34], and the improved time frequency representation algorithms are proposed in [35], which have been applied to mechanical fault signal analysis.

In view of the performance degradation of the conventional bispectrum methods in -stable distribution environment, the improved fractional low-order direct bispectrum and fractional low-order indirect bispectrum methods have been proposed for -stable distribution environment in this paper, and the improved fractional lower-order autoregression model bispectrum and fractional lower-order autoregressive moving average model bispectrum methods are presented for decreasing the estimate variance and increasing the bispectral flatness. We also summarize their calculation steps. The improved bispectrum methods and the traditional bispectrum methods employing second-order statistics are compared under Gaussian and -stable distribution environments; the simulation results show that the improved bispectrum methods have performance advantages compared to the traditional methods. Finally, we apply the improved methods to estimate the bispectrum of the normal and outer race fault signal; the result indicates that the proposed methods are feasible and effective for fault diagnosis.

In this paper, several improved bispectrum analysis methods based on fractional lower-order statistics are proposed for mechanical bearing fault diagnosis in Gaussian or -stable distribution noise environment. The paper is structured in the following manner. -stable distribution and the bearing fault signals are introduced in Section 2. The improved fractional lower-order bispectrum methods are demonstrated, and the simulation comparisons employing the traditional bispectrum methods and the improved bispectrum methods are performed to show the advantage of the proposed methods in Section 3. The conjoint application simulations of the actual bearing fault signals employing the proposed bispectrum and time-frequency distribution in [35] are demonstrated in Section 4. Finally, the conclusions and future research are given in Section 5.

2. Bearing Fault Signals

The actual bearing fault signals data are obtained from the Case Western Reserve University (CWRU) bearing data center [36]. The experimental equipment adopts 6205-2RS JEM SKF type bearing, the outer race diameter is 20.472 inches, and the inner race and the ball diameter are 0.9843 inches and 0.3126 inches, respectively. The bearing outer race thickness is 0.5906 inches, motor load is 0 HP, and motor speed is 1797 rpm. The bearing faults of inner race, ball, and outer race are set, and the fault diameters are all 0.021 inches. The fault data are collected at 12,000 samples per second, and the outer race position relative to load zone is centered at 6:00. The normal signals are given in Figure 1(a), and the fault signals of inner race, ball, and outer race are shown in Figures 1(b)1(d), respectively. We can know that the waveform of the fault signals has a certain impulse.

In order to further verify the pulse characteristics of the bearing fault signals, we use -stable distribution statistical model to estimate the parameters of the inner race fault, ball fault, and outer race fault signals, and the results are given in Table 1 [37, 38]. As it can be seen, the characteristic index of the normal signals is equal to 2, which is Gaussian distribution. However, the characteristic index of the bearing fault signals is greater than 1 but smaller than 2, and it belongs to non-Gaussian -stable distribution ().

PDFs of the signals of inner race fault, ball fault, and outer race fault are shown in Figures 2(a)2(c), respectively. From the PDFs of normal and fault signals, we can see that PDFs of fault signals have heavy tails. Most of the parameters are approximately equal to zero in Table 1, and Figure 2 shows that PDFs of the fault signals are near symmetric. Hence, distribution is a more concise and accurate statistical model for the bearing fault signals.

3. Fractional Lower-Order Nonparametric Bispectrum Methods

3.1. Fractional Lower-Order Direct Bispectrum Method

is samples of the observation data; its discrete fractional lower-order Fourier transform is defined as where is a given real constant and . denotes order moment of , when is a real signal, , , and when is a complex signal, , where denotes conjugate operation.

Dividing the data into segments, each segment is points, and two adjacent segments overlap points; then and . According to equation (1), the discrete fractional low-order Fourier transform of the th segment can be written as

The corresponding discrete fractional low-order Fourier transform coefficient is defined as

Let the sampling frequency of the signal be and the total frequency sampling points be ; then . Combining equation (3), fractional low-order triple correlation of can be expressed aswhere , and we let , , and . The corresponding fractional low-order direct bispectrum of the th segment is given bywhere and . Averaging the bispectrum of those segments, fractional low-order direct bispectrum of can be gotten.

Fractional low-order direct bispectral estimation process is to first calculate the discrete fractional low-order Fourier transform coefficient and then compute its fractional low-order triple correlation and finally average fractional low-order triple correlation of all segments.

3.2. Fractional Lower-Order Indirect Bispectrum Method

By using the fractional lower-order moment, we define discrete fractional lower-order three-order cumulants (FLOTOC) of aswhere , , and . If are real, then the estimation of FLOTOC is given byand if are complex, thenwhere and .

From equations (7)–(9), The th segment of FLOTOC of the signal can be given bywhere , averaging fractional low-order third-order cumulant of the segments; we can get

Taking the windowed two-dimensional discrete Fourier transform of equation (11), fractional low-order indirect bispectrum estimation of the signal can be given bywhere is a two-dimensional window function and . The bispectrum estimation process of fractional low-order indirect method is to first calculate the discrete fractional low-order third-order cumulant of each segment, followed by averaging the fractional low-order third-order cumulant of all segments, and finally compute the two-dimensional Fourier transform.

3.3. Application Review

In this simulation, the test signal is defined aswhere is three cosinoidal signals and is additive Gaussian noise or distribution noise. , , , , , , , and . When is additive Gaussian noise, signal to noise ratio () can be used. But is distribution noise, is inapplicable, and generalized signal to noise ratio () is written aswhere is the dispersion coefficient of distribution noise. According to the given , the amplitude of the signal is written as

Letting and , the traditional bispectrum direct and indirect methods and fractional lower-order bispectrum direct and indirect methods are applied to estimate the bispectrum of the signal under Gaussian distribution noise and distribution noise; the simulation results are shown in Figures 36.

3.4. Remarks

The direct bispectrum estimations of the signal x(n) under Gaussian noise environment () are shown in Figure 3. Figures 3(a) and 3(b) are the traditional direct bispectral estimation and its three-dimensional graph estimation, respectively. Fractional lower-order direct bispectral estimation and its three-dimensional graph are given in Figures 3(c) and 3(d), respectively. The result shows that both methods can estimate the bispectrum of the signal x(n) well. Figure 4 shows the direct bispectral estimations of the signal x(n) under distribution noise environment (; ). We can know that the traditional direct bispectral estimation fails in Figures 4(a) and 4(b), but the proposed fractional lower-order direct bispectral estimation method in Figures 4(c) and 4(d) shows good performance. Figures 5 and 6 are the traditional indirect bispectrum and fractional lower-order indirect bispectrum of the signal x(n) under Gaussian noise environment () and distribution noise environment (; ), respectively. The results show that both methods have better performance under Gaussian noise environment, but the conventional indirect bispectrum method degenerates under noise environment, and fractional lower-order indirect bispectrum method can better estimate out the bispectrum of the signal x(n). Hence, the fractional lower-order direct and indirect bispectrum methods are robust.

The fractional lower-order direct and indirect bispectrum methods have larger variance; we can increase the number of segments and segment length by overlapping adjacent segments to reduce variance and add two-dimensional window function to improve the frequency resolution of bispectral estimation. The proposed fractional low-order bispectral (FLOB) estimation methods still have biperiodicity and symmetry, so they can be used to quickly calculate the bispectrum of the signal.

Fractional lower-order bispectral biperiodicity properties are as follows:

Fractional lower-order bispectral symmetry properties are as follows:

4. Fractional Lower-Order Parametric Bispectrum Methods

4.1. Fractional Lower-Order AR Model Bispectrum Method
4.1.1. Principle

The conventional real order AR model process can be expressed as

The traditional AR model bispectrum of the signal in (18) is given bywhere is the estimate of the third moment of the driving noise, is the system transfer function, and

A fractional lower-order AR model process may be written aswhere is order of the AR model and are its parameters, which are real numbers. is an independent identically distributed (i.i.d) random process, is its characteristic index, and is its dispersion coefficient. can be expressed by a finite impulse response (FIR) system [27].where is system impulse response coefficient and is i.i.d random process. Taking z transformation for equations (21) and (22) and obtaining the system transfer function,

Taking -order moment with equation (23), we have, and on the unit circle; then

According to the definition of the AR model bispectrum in (19) and -order moment in (24) and (25), we define fractional lower-order AR model bispectrum (FLOARB) asand FLOARB on the unit circle is written aswhere is the dispersion coefficient of the driving random process and is conjugate operation.

According to the definition of the fractional lower-order three-order cumulants in (8), we define , , andwhere is named as fractional low-order three-order cumulants matrix (FLOTOCM), which is Toeplitz. From solving the AR bispectrum coefficients equations and fractional low-order moment matrix equations in [30], we let

Then,

Equation (30) is a fractional low-order bispectrum Yule-Walker equation. The AR model parameters can be gotten by solving equation (30). When are substituted into equation (27), we can get fractional low-order AR model bispectrum estimation.

In this paper, we apply the final prediction error (FPE) criteria to determine the order of fractional low-order AR model. When increases gradually from 1, FPE will be the minimum at a certain , which just is the most appropriate order. The calculation formula can be written aswhere is the variance of residuals.

We summarize the steps of the FLOAR model bispectrum method as follows:Step 1: computing fractional lower-order three-order cumulants of the signal with equations (7)–(9) and constructing fractional low-order three-order cumulants matrix in (28).Step 2: solving fractional low-order bispectrum Yule-Walker equation in equation (30) and getting the coefficients .Step 3: computing fractional low-order AR model bispectrum of by substituting into equation (19).

4.1.2. Application Review

In this simulation, the experimental signal is in equation (13). The performances of the traditional AR model bispectrum method and the improved FLOAR model bispectrum method are compared under Gaussian distribution noise () and distribution noise (; ), and the simulation results are shown in Figures 710.

In order to further verify the advantages of FLOAR model bispectrum method, we conduct comparative experiment on two methods under different when GSNR = 20 dB, and their parameter estimations are shown in Figure 11. When and GSNR changes from 14 dB to 24 dB, we compare the changes of the errors power with the AR and FLO-AR model bispectrum methods under -stable distribution noise environment, and the simulation is given in Figure 10.

4.1.3. Remarks

Figure 7 is the AR model and FLOAR bispectrum estimations of the signal x(n) under Gaussian noise environment (). The traditional AR model bispectral estimation and its three-dimensional graph estimation are shown in Figures 7(a) and 7(b), respectively. Fractional lower-order AR model bispectral estimation and its three-dimensional graph are given in Figures 7(c) and 7(d), respectively. The result shows that both methods can estimate the bispectrum of the signal x(n) well under Gaussian noise environment.

The AR and FLOAR model bispectrum estimations of the signal x(n) under distribution noise (; ) are given in Figure 8. The simulation result shows that the conventional AR bispectrum method fails under noise environment in Figures 8(a) and 8(b), but the improved FLOAR bispectrum method shows good toughness in Figures 8(c) and 8(d). As a result, the AR bispectrum method is only suitable to analyze the signals in Gaussian environment, but the FLOAR bispectrum method can work in Gaussian and noise environment, which is robust.

4.2. Fractional Lower-Order ARMA Model Bispectrum Method
4.2.1. Principle

A fractional lower-order autoregressive moving average (FLOARMA) model i.i.d process can be given bywhere and are orders of the AR and MA model, respectively. and are their parameters, which are real numbers. is an independent identically distributed (i.i.d) random process, is its characteristic index, and its dispersion coefficient is .

Simplifying equation (32), we havewhere and . Taking z transform on both sides of equation (33), we obtainwhere is system impulse response coefficient, is system transfer function, and is i.i.d random process. The AR model parameters changes with ; and are the FIR and IIR filters, respectively.

Taking -order moment with equation (36),, and on the unit circle; then

According to the definition of the ARMA model bispectrum and -order moment of in (38) and (39), we define fractional lower-order ARMA model bispectrum (FLOARMAB) asand FLOARB on the unit circle is written aswhere is the dispersion coefficient of the driving random process and is conjugate operation.

To obtain the coefficients of the fractional lower-order ARMA model and , we should multiply on both sides of equation (32) taking fractional lower-order covariance to getwhere

It can be shown [24] that, for the fractional lower-order covariance of the signal and noise in (45),

Then, if , we have

Equations (47) and (48) are a generalized Yule-Walker equation, and we can obtain the fractional lower-order AR model (FLOARM) coefficients by solving them.

If , equation (42) changes as

Letting , we have

We can obtain the coefficients of the fractional lower-order MA model by solving the nonlinear equation (35).

A finite -order FLOMA model can be equivalent by an approximate infinite -order FLOAR model; we havewhere is error. Letting , we obtain

We should multiply on both sides of equation (52). We take fractional lower-order covariance to getwhere and . We can also obtain the coefficients of the FLOMA model by solving the Yule-Walker equation (54).

We summarize the steps of the FLOARMA model bispectrum method as follows:Step 1: solving fractional low-order bispectrum Yule-Walker equation in (48) and getting the coefficients of the FLOAR model.Step 2: obtaining the coefficients of the FLOMA mode by solving the nonlinear equation (50) or the Yule-Walker equation in (54).Step 3: determining the order and employing FPE criterion.Step 4: computing fractional low-order ARMA model bispectrum of by substituting and into equation (40) or equation (41).

4.2.2. Application Review

In this simulation, in equation (14) is used as the experimental signal. The ARMA model and FLOARMA model bispectrum methods are compared to demonstrate their performance under Gaussian distribution noise () and distribution noise (; ), and the simulations are shown in Figures 912.

4.2.3. Remarks

The ARMA model and FLOARMA model bispectrum estimations of the signal x(n) are shown in Figure 9 under Gaussian noise environment (). Figures 9(a) and 9(b) are the existing ARMA model bispectral estimation and its three-dimensional graph estimation, respectively. FLOARMA model bispectral estimation and its three-dimensional graph are given in Figures 9(c) and 9(d), respectively. The results show that both methods can estimate out the bispectrum of the signal x(n) well under Gaussian noise.

Figure 10 shows the ARMA and FLOARMA model bispectrum estimations of the signal x(n) under distribution noise (; ). The simulation results show that the existing ARMA bispectrum method in Figures 10(a) and 10(b) fails under noise environment, but the improved FLOARMA model bispectrum method in Figures 10(c) and 10(d) has good toughness. Hence, the ARMA bispectrum method only works for Gaussian noise environment, but the FLOARMA bispectrum method can be applied in Gaussian and noise environment, which is robust.

5. Application Simulations

The impulse of the outer race fault signals in the vibration position of DE, FE, and BA is generated because of the local defects of rolling element bearings, and the waveforms are given in Figure 3(d) and Table 1. We can know that the fault signals are nonstationary and non-Gaussian -stable distribution process. In this section, the experiment signal adopts from the normal signal and the bearing outer race fault signal in the vibration position of DE, and 0.2-second data is selected as the test signal, which is collected at 12,000 samples per second; then . We apply the improved FLODB, FLOIDB, FLOAR, and FLOARMA bispectrum methods to analyze the normal and bearing outer race fault signals, and the simulations are shown in Figures 11 and 12. In this section, we have only extracted the first quadrant of the bispectral representation (, ) to analyze the signals.

Figures 11(a), 11(c), 11(e), and 11(g) are the bispectral contour maps of the bearing normal signal employing the FLODB, FLOIDB, FLOAR, and FLOARMA bispectrum methods, respectively. Their bispectral three-dimensional diagrams are shown in Figures 11(b), 11(d), 11(f), and 11(h), respectively. It is observed that the improved methods can effectively suppress the noise and estimate out the frequency components of the normal signal, and the spectral peaks exist near the central frequency of 1060 Hz; hence, the transient harmonic vibration components of the normal signal are about 1060 Hz.

Figures 12(a), 12(c), 12(e), and 12(g) are the bispectral contour maps of the bearing out race signal employing the FLODB, FLOIDB, FLOAR, and FLOARMA bispectrum methods, respectively, and their bispectral three-dimensional diagrams are shown in Figures 12(b), 12(d), 12(f), and 12(h), respectively. It is observed that the spectral peaks of Figures 12(a), 12(c), and 12(e) exist near the central frequencies of 600 Hz and 2800 Hz, and those in Figure 12(g) exist near the central frequencies of 600 Hz, 2800 Hz, and 3500 Hz. Hence, the outer race signal has fault harmonic vibration components of 600 Hz, 2800 Hz, and 3500 Hz, indicating that the bearing outer race is damaged. Figure 12(i) clearly shows the gap between the impacts regularly changing. The interval between impulses A, B, C, D, E, and F is approximately 30 ms; then the characteristic frequency of the bearing outer race fault is about 33.333 Hz. Figure 12(j) shows wave shape and spectrum of the envelope signal of the outer race fault signal using resonance-demodulation approach; we can clearly see the pulse separation. We can see that FLODB and FLOIDB methods have larger variance; the FLOARMA bispectrum method has lower variance and higher frequency resolution, and its performance is optimal.

To further verify the advantages of the proposed fractional low-order bispectrum methods, distribution noise (; ) is added in the -stable distribution outer race fault signal as the actual working environment background noise. The conventional bispectrum methods including the direct bispectrum, indirect bispectrum method, AR model, ARMA bispectrum methods, and the improved fractional lower-order bispectrum methods including FLODB, FLOIDB, FLOARB, and FLOARMAB methods are used to analyze the outer race bearing fault signal under distribution background noise environment. The simulations are given in Figures 13 and 14. The simulation results show that the existing bispectrum methods in Figures 13(a)14(h) fail, but the proposed fractional low-order bispectrum methods in Figures 14(a)14(h) have good performance, and it is observed that the transient harmonic vibration components are near the central frequencies of 600 Hz and 2800 Hz.

6. Conclusions

The bearing fault signals are a non-Gaussian and nonstationary process, and -stable distribution is a more appropriate statistical model for them. The improved FLODB, FLOIDB, FLOAR, and FLOARMA model bispectrum methods have been proposed for the fault signals employing fractional low-order statistics. The improved methods can be applicable to Gaussian and -stable distribution noise environment, and their performances are superior to the existing direct bispectrum method, indirect bispectrum method, and AR and ARMA model bispectrum analysis methods. Fractional low-order nonparametric bispectrum estimation methods, FLODB and FLOIDB, require a large number of data samples and have a large estimation variance, but the fractional low-order parametric bispectrum estimation, FLOAR and FLOARMA model bispectrum, has small variance and produces fewer parameters that describe the characteristics of the target; hence, it can be directly used for target features. We can apply the improved methods to analyze the -stable distribution bearing fault signal, even -stable distribution noise environment, and the fault characteristic frequency, the dominant frequency, and the other fault frequency features of the fault signals can be clearly obtained. Combining the fractional low-order time-frequency methods, more fault characteristics can be obtained, and the joint diagnosis will be realized for the bearing fault signals. In the future, we can also apply the bispectrum diagonal slice to reflect the coupling information between fault signals so as to realize the recognition of the fault characteristics. The complete mechanical bearing fault state spectrum can be established based on fractional low-order bispectrum estimation for the fault signals, which can provide a new way for the fault diagnosis and online monitoring.

Data Availability

The supplementary file in txt file format is the original experimental data of the paper.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was financially supported by the National Natural Science Foundation of China (61261046), Natural Science Foundation of Jiangxi Province China (20192BAB207002), Science and Technology Project of Provincial Education Department of Jiangxi (GJJ170954), Science and Technology Project of Jiujiang University, China (2014SKYB009), and Science and Technology Project of Jiangxi Provincial Health Commission (SKJP220200278).

Supplementary Materials

(1) The txt file ball fault signals i the signal when the ball fails. The ball diameter is 0.3126 inches, and the fault data are collected at 12,000 samples per second. (2) The txt file inner race fault signal is the signal when inner race fails. The inner race diameter is 0.9843 inches, and the fault data are collected at 12,000 samples per second. (3) The txt file normal signals are a trouble-free signal, and the fault data are collected at 12,000 samples per second. (4) The txt file outer race fault signals is the signal when the outer race fails. The bearing outer race thickness is 0.5906 inches, and the fault data are collected at 12,000 samples per second. (Supplementary Materials)