Bearing incipient fault diagnosis based upon maximal spectral kurtosis TQWT and group sparsity total variation denoising approach

Abstract. Localized faults in rolling bearing tend to result in periodic shocks and thus arouse periodic responses in the vibration signal. In this paper, a novel fault diagnosis method based on maximal spectral kurtosis tunable Q-factor wavelet transformation (TQWT) and group sparsity total variation denoising (GS-TVD) is proposed to address the issue of bearing incipient failure. Firstly, the range of Q-factor was pre-selected according to the spectral distribution of impulse component, and bearing vibration signal was transformed by the TQWT method. Then, the spectral kurtosis of each scale transform coefficients was calculated, and the optimal Q-factor and decomposition scale can be selected according to the kurtosis maximum principle. In order to remove the interference components and high-frequency noise from the reconstructed vibration signal generated by inverse TQWT, the GS-TVD approach is employed, thus the cyclic periodicity characteristic and transient impulses can be detected obviously. The two cases experimental results indicate that the proposed technique is more effective and applicable for bearing incipient fault diagnosis compared with traditional method.


Introduction
Rolling element bearings have been widely applied in rotating machines of the modern industry.During operation, some common and unknown failures such as pitting, spalling and localized wear are inevitable, and these bearing faults seriously affect the safety of machinery.Therefore, bearing fault detection and diagnosis have attracted substantial attention in recent years [1][2][3].
Vibration signal analysis is most widely used to diagnose rolling bearing faults.When a bearing operates with localized defect, it would excite the resonance frequency modulated by characteristic frequencies caused by bearing defects.This would cause the measured vibration signal with the nonlinear and non-stationary characteristics [4].However, at the early stage, the fault symptoms generated by a localized damage has low amplitude, which are often immersed in heavy background noises and rotation frequencies.Common used indexes such as root mean square (RMS), peak-to-peak value and Kurtosis value, are effective to describe the degradation process of a bearing.However, the bearing fault type cannot be identified easily using these indexes and the certainty of the fault occurrence requires to be double assured.
In the past years, many vibration signal processing techniques, such as wavelet packet transform [5], ensemble empirical mode decomposition (EEMD) [6], local mean decomposition (LMD) algorithm [7] and Variational mode decomposition (VMD) algorithm [8], higher order energy operator fusion [9], Quaternion singular spectrum [10], blind source separation method [11], low-rank matrix approximations [12], sparse representation [13,14] and deep learning [15], etc., have been developed to extract the fault features from measured signals for bearing fault diagnosis.Currently, the turntable -factor wavelet transform (TQWT) was original proposed by Selesnick, and the advantage of TQWT is that the -factor is easily and continuously adjustable [16,17].When a fault occurs on bearings, the measured vibration signals contain both the fault transient impact component and sustained oscillatory component.Therefore, the TQWT has the high -factor and low -factor simultaneously, which is very suitable to separate the fault related components from the sustained oscillatory components.Ref. [18] recently proposed ensemble empirical mode decomposition (EEMD) and TQWT method to extract the weak fault features of rolling bearing in an accelerated bearing life test.In this method, the EEMD is applied to decompose the vibration signals into a series of IMFs and then the TQWT is employed to separate the main IMF into high -factor component and low -factor component to diagnose the early fault.Furthermore, an improvement works has been aroused in applying intrinsic characteristic-scale decomposition (ICD) and TQWT for fault diagnosis of rolling bearings [19].However, there are two main shortcomings in the above methods: 1) The EEMD and ICD methods are time consuming in signal decomposition, and also the self-adaptive, mode mixing phenomenon and end effect problems cannot be eradicated.
2) Although the TQWT is an interesting attempt to extract the fault characteristics, the selection of the proper -factors is still a problem.As the filtered signal processed by -factors is far away from its theoretical central frequency, and the bandwidth is often heavily polluted, it is difficult to represent the factual features clearly.
To address the above issues, in this paper, a novel fault diagnosis method based on group sparsity total variation denoising (GS-TVD) and maximal spectral kurtosis TQWT is proposed.The optimal -factors and decomposition level can be found for fault feature extraction by the maximal spectral kurtosis TQWT, which match the fault periodical impulses well.Unfortunately, it should be noted that some interference components and high-frequency noise may still be present in the optimal -factors and decomposition level area.To overcome this problem, unlike the traditional signal decomposition method, the group sparsity total variation denoising (GS-TVD) approach is developed for noise reduction, which can preserve sharp edges in the underlying signal so that the transient impulse will be preserved for feature extraction.Finally, results of incipient feature extraction indicate that the period of transients can be detected more accurately and effectively in the case that previous approach failed, which can significantly improve the performance of proposed method for extracting transient impulses from heavy noisy vibration signal.
The rest of this paper is organized as follows.Section 2 mainly introduces the maximal spectral kurtosis TQWT method and group sparsity total variation denoising method in detail.The flow chart of the proposed model for bearing incipient fault diagnosis is presented in Section 3. Section 4 investigates the effectiveness of the proposed method compared with the other method in feature extraction by using two bearing fault cases.Conclusions are drawn in Section 5.

Maximal spectral kurtosis tunable 𝑸-factor wavelet transform
The resonance characteristics of oscillatory signal can be depicted by quality factor  ,  =   ⁄ , where  represents center frequency of signal and  is bandwidth.The tunable -factor wavelet transform (TQWT) was proposed by Selesnick as a flexible discrete wavelet transform for oscillatory signal processing [16].The TQWT method includes three easily changeable parameters, i.e. quality factor , total over-sampling rate  and number of levels of decomposition .For each level of TQWT decomposition, the signal () with sampling frequency  can be decomposed into low-pass and high-pass sub-band signals with sampling frequencies  and , respectively, the parameters  and  are scaling values.Moreover, the low-pass filter  () and low-pass scaling  are used to obtain low-pass sub-band.Similarly, the high-pass sub-band can be generated by high-pass filter  () and high-pass scaling  accordingly.The TQWT approach uses following given low-pass and high pass filters: where () is the frequency response of the Daubechies filter which has two vanishing moments.
The () can be defined as follows: The -factor  and parameters  can be formulated via filter bank parameters  and  as follows [16]: where  represents center frequency of signal and  is bandwidth of sub-band .For more details about decomposition and reconstructed operations can be found from Selesnick [16].
For maximal spectral kurtosis TQWT method, it applies the principle of maximum kurtosis, the kurtosis of each scale wavelet coefficient can be calculated once the original vibration signals are decomposed by TWQT method, then the calculated Kurtosis value can be used for selecting the optimal frequency band, finally, the corresponding reconstructed signal is obtained by inverse TQWT method.

Group sparsity total variation denoising
Total variation denoising (TVD) assumes that the noisy signal  is of the form [20][21][22]: where  is unknown signal and  white Gaussian noise.The unknown signal  can be obtained by solving the following optimization problem: The regularization parameter  > 0 controls how much smoothing is performed, and matrix  is first-order difference matrix, i.e.: Note that  =  and penalty function () can be defined by: Here, the parameter  represents the group size.
If  = 1, the penalty function () = ‖‖ and Eq. ( 6) is the standard 1-D TVD problem.If  > 1, the penalty function () is a convex measure of group sparsity.In this paper, we discuss group-sparse total variation denoising (GS-TVD) method.
Definition: A -point group can be denoted by the vector : Thus, the penalty function () can be modified and written by: To calculate the optimization of (), the optimization of penalty function () is calculated firstly.Not that for all  and  ≠ 0 with equality when  = , we have: The Eq. ( 9) is used for each group, thus the optimization of () is presented as follows: with (, ) ≥ (), (, ) = ().Note that (, ) is quadratic equation on .The Eq. ( 10) can be written as: where  does not depend on the value of , and Λ() is a diagonal matrix, i.e., Λ() = diag(| |).
In order to minimize the function (), the majorization-minimization (MM) algorithm [23] defines an iterative algorithm via: where  is iteration exponent.The iteration is initialized with  ( ) .Thus the Eq. ( 14) is written as: Solving the above equation, we have: Note that the problem with update Eq. ( 16) is that as the iterations progress, some values of Λ( ( ) ) will generally go to zero and becomes sparse, and therefore some entries of Λ ( ( ) ) in Eq. ( 16) will go to infinity, leading the result is not accurate.We use the matrix inverse lemma from Ref. [24] to address this issue.Generally, the matrix inverse lemma has several forms, a common form is: Thus, the [ +  Λ( ( ) )] in Eq. ( 16) can be written as: Finally, using Eq. ( 18), the update Eq. ( 16) can be written as: Therefore, the group-sparse total variation denoising (GS-TVD) problem can be solved by the iterative algorithm in Eq. (19).
In order to investigate and compare the denoising performance of GS-TVD and general TVD method, Fig. 3 illustrates group-sparse TV denoising on a blocks signal.The test block signal and its corrupted signal are provided in Fig. 3(a).Fig. 3(b) and (c) respectively illustrate the denoising performance of GS-TVD and general TVD method.Moreover, denoising errors are shown in Fig. 3(d).By comparison, the simulation results indicate that the significant improvement in denosing is obtained with GS-TVD approach compared to the TVD method.Here, the group size was set to  = 2, the regularlization parameter is 2, and the number of iteration is 20.
a) The test signal and noisy test signal b) The test signal and denoising signal using TVD method c) The test signal and denoising signal using GS-TVD method d) Denoising error using TVD and GS-TVD method Fig. 3. Comparison of simulation results between GS-TVD compared to the TVD method

Feature extraction of bearing fault using maximal spectral kurtosis TQWT and GS-TVD approach
In this paper, a new incipient fault feature extraction technique by using maximal spectral kurtosis TQWT and group sparsity total variation denoising is proposed to extract fault frequencies from its measured signal.The procedures of fault feature extraction are divided into 2) Calculate the kurtosis of each scale wavelet coefficient obtained in step (1) and select the local maximum area according to its maximum kurtosis value.Keep the scale coefficient of the maximum area unchanged, and the other scale coefficients are set to 0.
3) Apply the inverse TQWT method to the scale coefficient of the maximum area and the corresponding reconstructed signal is obtained.
4) Determine the appropriate group size, and utilize the group sparsity total variation denoising method to reduce noises in the reconstructed signal, and the corresponding filtered signal is obtained.
5) Calculate envelope spectrum of the filtered signal.6) Determine fault types by contrasting envelope spectral peak and the theoretical fault character frequency.
The flow chart of the proposed fault feature extraction technique is illustrated in Fig. 4.  1.
Fig. 6 shows the Kurtosis curve over the whole life-cycle of bearing 1 and shows that there is a long time in stable or normal operation for the bearing in whole life-cycle and the period of fault occurrence and severity is relatively short.As shown in Fig. 6, there is an obvious transient feature for the incipient fault at the stage of point 647.However, due to the interference of environment and background noises, the engineers are not sure whether the fault is happened before point 647 or not.Hence, to verify the effectiveness of the proposed method for bearing incipient fault diagnosis, the experimental data at point 535 was chosen which has no obviously wave phenomenon during the whole life-cycle.

Results and discussion for Case 1
The original vibration signal (2048 sampling points is selected, namely, approximately 0.1 s), the amplitude spectrum and its corresponding envelope spectrum of outer-race at point 535 are displayed in Fig. 7 Fig.7(c), although the spectrum peak at 230 Hz which consists with the outer-race fault frequency can be detected without denosing, however, the spectrum peak masked by heavy background noise and features are not be evident enough to detect fault.The proposed method is then employed to analyze the outer-race vibration signal.First, for parameter  setting, it can be seen from the amplitude spectrum (see Fig. 7(b)) that the spectrum of the shock component is wider, and if the value of  is selected higher, due to that the noise component will increase and the underlying wavelets will have more oscillations with a narrower frequency responses, hence, we set the scope of resonant quality factor is  ∈ [1,2] with the interval ∆ = 0.1.The maximal spectral kurtosis TQWT method is applied on the original vibration signal, then the kurtosis values under different scales- are obtained.The Kurtosis surf plot under different -factor and decompose level is displayed in Fig. 8(a).It can be seen from Fig. 8(a), the maximum index of kurtosis 4.709 is emerged in -factor  = 1.55 and level  = 3.5, as marked by the black rectangle, which means that this area contains much more transient impulse components information than other  -factor and levels areas.Applying the inverse TQWT method to this area and the corresponding reconstructed fault signal is shown in Fig. 8(b).Obviously, the reconstructed signal illustrated in Fig. 8(b) shows that the background noise is suppressed relatively.However, the low-frequency noise components still exist, mainly because of the range of frequencies available band.By computation, the minimum frequency is: and the maximum frequency is  =   /2 ≈ 4688 Hz, thus the impact component with frequency 4252 Hz (see Fig. Furthermore, the group-sparse TVD method is employed to process the reconstructed bearing vibration signal, here, the group size is set to  = 3, the regularization parameter  = 2 and the iteration time is 300.The filtered signal with group-sparse TVD method and its envelope spectrum are illustrated in Fig. 9 in reconstructed vibration signal have been removed evidently as shown in Fig. 9(a), and the cyclic periodicity characteristic is prominent in the filtered vibration signal.From Fig. 9(b), an obvious spectrum line 230.1 Hz, which is close to the  frequency according to the Table 1, and its harmonics are exposed in the envelope spectrum, meanwhile, the amplitudes of fault frequencies illustrated in Fig. 9(b) are evidently enhanced compared with that in Fig. 8(c), thus making the feature extraction becomes easier and accurate.Moreover, besides the outer race fault frequency, it is observed that the spectrum peaks at 70 Hz and 100 Hz are also identified in Fig. 9(b), which consists with the double shaft rotating frequency (70 Hz ≈ 2  = 33.33×2)and its triple frequencies (100 Hz ≈ 3 = 33.33×3).The above results illustrate that the fault characteristic frequencies can be detected by the proposed method even if the useful failure information is extremely weak.To verify the superiority of the proposed approach, the Spectral Kurtogram (SK) combined with EEMD method is employed to analyze the original vibration signal.The SK is a powerful analysis tool to detect transients feature from their noisy vibration signals and has been widely studied and applied in the rotating machine diagnosis, see refs [27] and [28].First, we apply the EEMD method [29,30] to decompose the original vibration signal, thus the 12 IMF components are generated and displayed in Fig. 10, which can be helping to distinguish the periodic impulse from the mixed noisy signal.Then, we calculate the Kurtosis values of all the IMFs model orderly, see Table 2, and 2nd IMF model that contains the most fault information according to its maximum kurtosis value.The Kurtogram of 2nd IMF model is displayed in Fig. 11(a), from which the optimal demodulation frequency band namely [0 Hz-2560 Hz] can be detected, as marked by the black rectangle.Thus, a band-pass filter is designed to extract the potential features from the 2nd IMF model.Last, the envelope spectrum is applied to the filtered signal, the corresponding envelope spectrum is illustrated in Fig. 11(b).As can be seen, the fault characteristic frequencies  = 236.4Hz cannot be detected in the envelope spectrum and it is also hard to distinguish the fault location from the incipient vibration signal.In conclusion, the comparison above re-confirms the efficiency of the proposed approach.

Case 2 -experimental setup
The bearing vibration signals were generated by the Curtin University in Australia [31].The

Results and discussion for case 2
In this experiment, the inner race vibration fault signal (5120 points), i.e. 0.1 second data was selected.In order to get a real early failure, a Gaussian white noise with amplitude  = 17 (a high amplitude noise was added for highlighting the effect of proposed method) was added to the original signals.The time domain waveform, amplitude spectrum and its corresponding envelope spectrum are shown in Fig. 13(a), (b) and (c), respectively.However, it can be found that the periodical pulses are submerged in heavy noise in the time domain signal, and the inner fault  cannot be identified from the envelope spectrum and fault type cannot be determined yet.
The proposed approach is applied to the experimental signal, the algorithm parameters of Similarly, the reconstructed bearing vibration signal is also processed by the group-sparse TVD method, and the group size is set to  = 2, the regularlization parameter is  = 2 and the iteration time is 300.The filtered signal with GS-TVD method and its wavelet time-frequency diagram, and its envelope spectrum are displayed in Fig. 15.
Obviously, it can be observed that the GS-TVD algorithm not only extracts the transient impulse components but also the noise components in reconstructed vibration signal have been removed evidently.As shown in Fig. 15(a), the phenomenon of frequency modulation can be clearly observed by the wavelet time-frequency diagram, and the impact interval of 9.7 ms could be also identified via the wavelet time-frequency diagram and the time-domain waveform, which coincides with the inner race fault frequency.From Fig. 15(b), as can be seen, the characteristic frequency 100.2 Hz, which is agreeable with the theoretical fault frequency of inner race 103.59 Hz, and its harmonic frequencies are clearly detected, hence, the proposed approach has a very powerful effect on the incipient fault diagnosis.shown in Fig. 16.Meanwhile, the Kurtogram is also calculated to assist selecting the optimal frequency band for envelope demodulation analysis, as illustrated by Fig. 17(a).Based on highest SK value in Kurtogram, see Table 4, an optimal frequency band namely [0-1066.66Hz] is selected for the 4th IMF model, as marked by the black rectangle, then the sensitive component within potential fault frequency is separated via a band pass filter.Unfortunately, the inner race fault characteristic frequencies  = 103.59Hz cannot be detected in the envelope spectrum, see Fig.   Obviously, although the double inner fault frequency (210 Hz) can be identified, the peak amplitude is lower than that of the proposed approach.Therefore, the comparison analysis demonstrates that the proposed maximal spectral kurtosis TQWT and group sparsity TVD approach can not only provide better de-nosing performance but also has superiority in enhancing performance for incipient fault diagnosis.

Conclusions
In this paper, a novel maximal spectral kurtosis TQWT and group sparsity total variation denoising (GS-TVD) technique is proposed for bearing incipient fault diagnosis.Two main contribution of this work can be summarized as: 1) The optimal -factors and decomposition level can be found for fault feature extraction by the maximal spectral kurtosis TQWT, which match the fault periodical impulses well, overcoming the disadvantage of not being able to choose automatically using traditional TQWT method.
2) The group sparsity total variation denoising (GS-TVD) approach is applied to remove the noise from the reconstructed vibration signal generated by maximal spectral kurtosis TQWT, and the periodical transient impulses were detected obviously.Experimental results show that the proposed approach has better performance than existing algorithm.
Amplitude/m•s-2   Denoising error using TVD Denoising error using GS-TVD six steps:1) Preselect the scope of quality factor- according to the spectral distribution of impulse component, and the bearing vibration signals are decomposed by TWQT within this scope.

Fig. 4 . 4 . 1 .
Fig. 4. Flow chart of the proposed model for bearing incipient fault diagnosis 4. Experimental evaluation 4.1.Case 1 -experimental setup The bearing fault vibration signals were generated by the NSFI/UCR Center for Intelligent Maintenance Systems (IMS) [25, 26].The experimental setup is shown in Fig. 5.The rotation speed is 2000 rpm and the sampling rate is set at 20 kHz.Vibration acceleration signals were collected every 10 mins by a NI 6062E-DAQ-Card.The accelerated life test experiment was carried out successively for 8 days until the magnetic plug exceeds a certain level and causes an electrical switch to close.Meanwhile, the severe wear failure of outer race was detected in

Fig. 5 .Fig. 6 .
Fig. 6.The Kurtosis curve of whole life-cycle of the rolling bearing 1 (a), (b) and (c), respectively.From Fig. 7(a), the cyclic periodicity characteristic of time domain waveform relevant to bearing fault cannot be discovered.Moreover, as shown in 2863.BEARING INCIPIENT FAULT DIAGNOSIS BASED UPON MAXIMAL SPECTRAL KURTOSIS TQWT AND GROUP SPARSITY TOTAL VARIATION DENOISING APPROACH.QING LI, STEVEN Y. LIANG

Fig. 7 .
7(b)) and high-frequency noise were contained within the scope [1011 Hz, 4688 Hz].This conclusion can be verified by the corresponding envelope spectrum of the reconstructed signal as illustrated in Fig. 8(c).From Fig. 8(c), although the outer race fault frequency ( ) can be identified, several unrelated frequencies around  still remain in frequency domain and also the  information are not very obvious.a) Original experimental signal b) Amplitude spectrum of the original bearing vibration signal c) Envelope spectrum of the original bearing vibration signal Original experimental signal, amplitude spectrum and envelope spectrum (a) and Fig. 9(b), respectively.It can be observed that the noise components 0 INCIPIENT FAULT DIAGNOSIS BASED UPON MAXIMAL SPECTRAL KURTOSIS TQWT AND GROUP SPARSITY TOTAL VARIATION DENOISING APPROACH.QING LI, STEVEN Y. LIANG 1418 JOURNAL OF VIBROENGINEERING.MAY 2018, VOLUME 20, ISSUE 3

Fig. 8 . 9 .
a) Kurtosis surf plot under different -factor and decompose level b) The reconstructed vibration signal c) The envelope spectrum of the reconstructed vibration signal Kurtosis surf plot under different -factor and decompose level, and the reconstructed vibration signal from the maximum peak a) The filtered signal with GS-TVD method b) The envelope spectrum of the filtered signal Fig.The filtered bearing signal and its envelope spectrum processed by proposed method INCIPIENT FAULT DIAGNOSIS BASED UPON MAXIMAL SPECTRAL KURTOSIS TQWT AND GROUP SPARSITY TOTAL VARIATION DENOISING APPROACH.QING LI, STEVEN Y. LIANG a) IMF1-IMF6 b) The Hilbert envelop spectrum of IMF1-IMF6 c) IMF7-IMF12 d) The Hilbert envelop spectrum of IMF7-IMF12 Fig. 10.Sub-band components of original signal decomposed by using EEMD

11 .
INCIPIENT FAULT DIAGNOSIS BASED UPON MAXIMAL SPECTRAL KURTOSIS TQWT AND GROUP SPARSITY TOTAL VARIATION DENOISING APPROACH.QING LI, STEVEN Y. LIANG experimental setup (Machinery Fault Simulator, MFS) is shown in Fig. 12.The rotation speed is 29 Hz, sampling frequency is 51200 Hz and length of record is 10 seconds.The bearing information and the fault frequencies of inner race and outer race are summarized in Table 3.In this work, the inner race fault signal was chosen as an example for the subsequent analysis, and the theoretical inner race fault frequency is 103.59Hz based on the equation  = 1 + cos , where  and  are ball diameter and pitch diameter of the rolling bearing, respectively.Parameter  represents contact angle ( = 0),  denotes the number of elements,  represents shaft rotating speed ( = 29×60 = 1740 rpm).a) Kurtogram of second IMF model component b) The envelope spectrum of band-pass filtered signal Fig. Kurtogram of second IMF model component and its frequency spectrum

17 .
a) Kurtogram of 4th IMF model component b) The envelope spectrum of band-pass filtered signal Fig. Kurtogram and its frequency spectrum max =7 @ level 4.5, Bw= 1066.6667Hz,

Table 1 .
Geometric parameters and expected characteristic frequencies of the tested bearing

Table 3 .
Bearing information and fault frequencies (note that SS is shaft speed, i.e., 29 Hz) Bearing No. Number of balls Ball diameter Pitch diameter 13r a comparison, the EEMD result, i.e.,13IMF components related to the inner race fault are BEARING INCIPIENT FAULT DIAGNOSIS BASED UPON MAXIMAL SPECTRAL KURTOSIS TQWT AND GROUP SPARSITY TOTAL VARIATION DENOISING APPROACH.QING LI, STEVEN Y. LIANG