Intelligent Mechanical Fault Diagnosis Based on Multiwavelet Adaptive Threshold Denoising and MPSO

1 School of Mechanical Engineering, Jiangnan University, 1800 Li Hu Avenue, Wuxi, Jiangsu 214122, China 2 Jiangsu Key Laboratory of Advanced Food Manufacturing Equipment and Technology, Wuxi 214122, China 3 School of Mechanical & Electrical Engineering, Beijing University of Chemical Technology, Chaoyang District, Beijing 100029, China 4Graduate School of Bioresources, Mie University, Mie 514-8507, Japan


Introduction
Rolling element bearings are an important part of and widely used in rotating machinery.In practical application, bearing failures may cause the breakdown of equipment, and further, serious consequences may arise due to the failure.Thus, fault diagnosis and condition discrimination of bearings have an important significance for safe operation, guaranteeing production efficiency and reducing maintenance cost.Many reliability survey papers deal with failure statistics of rotating machinery subassemblies, focusing mainly on roller bearing because of their widespread use in industry [1][2][3][4].Occurrence rate of bearing faults is very high in rotating machines, and other faults arising in rotation machines are often associated with bearing faults.In many instances, the accuracy of the instruments and devices used to monitor and control the rotation machines is highly dependent on the dynamic performance of bearings.Although fault diagnosis of rolling bearings is often artificially carried out using time or frequency analysis of vibration signals, there is a need for a reliable, fast automated diagnosis method.
Vibration diagnosis is commonly used to detect the faults and identify the states in rotating machine.The condition diagnosis of rotating machinery depends largely on the feature analysis of vibration signals measured for the condition diagnosis because the signals carry dynamic information about the machine state [5][6][7].However, feature extraction for fault diagnosis is difficult, because if the vibration signals are measured at an early stage of the machine failure, or at a location away from the fault part, the vibration signals 2 Mathematical Problems in Engineering contain strong noise.Stronger noise than the actual failure signal may lead to misrecognition of useful information for the condition diagnosis.Thus, it is important that the feature of the signal can be sensitively extracted at the state change of a machine.
Wavelet transform (WT) is well known for its ability to focus on localized structures in time-frequency domain which has been widely used for fault diagnosis of rolling element bearings [8][9][10].It has the local characteristic of time domain as well as frequency domain and its time-frequency window is changeable.In the processing of nonstationary signals it presents better performance than the traditional Fourier analysis.However, the measured signals often contain strong noise and the fault features are hidden in the background noise, and it is not the best way for WT to match the different fault features with a single wavelet and scaling functions, which will reduce the fault diagnosis accuracy.Multiwavelet transform is the new development of WT.It is constructed from translations and dilations of scaling and wavelet vector functions and has the predominant properties such as orthogonality, symmetry, compact support, and higher order vanishing moments.Multiwavelet transform decomposes the signal into subsignals of different frequency bands based on vector basis functions, via inner product principle.Because of the multiple scaling and wavelet basis functions, multiwavelet transform has predominant advantages in feature extraction of signals.Recently multiwavelet transform has been applied in fault diagnosis of rotating machinery as a powerful tool.In [11], multiwavelet system was introduced to diagnose gear faults.In [12,13], multiwavelet lifting scheme was improved for compound faults separation and extraction.In [14], the undecimated multiwavelet was proposed for fault diagnosis of planetary gearboxes.
PSO algorithm is a population based stochastic optimization technique developed by Kennedy and Eberhart in 1995 and inspired by social behavior of bird flocking or fish schooling [15].In PSO algorithm, particles cooperate in finding good solutions for difficult discrete optimization problems.PSO algorithm has been applied to a variety of different problems, such as function optimization [16], scheduling [17], traveling salesman problem [18], neural network training [19,20], and clustering task [21][22][23] which is the topic of interest in this paper.In recent years, PSO algorithm has been successfully applied in mechanical fault diagnosis; the domestic research on PSO fault diagnosis issues also has many articles reporting [24][25][26][27].In [24], Bocaniala and Sa da Costa compared the time spent by PSO algorithm and genetic algorithm, testifying PSO algorithm with prominent superiority through fault diagnosis benchmark problem.In [25], Pan et al. used PSO algorithm to extract fault characteristics of rotation machinery.In [26,27], PSO algorithm was used to diagnose gearbox fault.In this study, a clustering model is constructed by using an improved PSO called MPSO algorithm.It is used to classify the SPs calculated from the signals in each machine state for condition diagnosis, as well as obtaining their optimal clustering centers.According to these optimal clustering centers' information, the conditions of the machine can be accurately identified.
In order to extract the fault features of signals more effectively and identify mechanical condition more accurately, this paper proposes a novel fault diagnosis method for rotation machinery based on multiwavelet adaptive threshold denoising and MPSO algorithm.GHM multiwavelet is employed for extracting weak fault features under heavy background noise, and the method of adaptively selecting appropriate threshold values for multiwavelet with energy ratio of multiwavelet coefficient is presented.The six nondimensional SPs in the frequency domain are defined to reflect the features of the vibration signals measured in each state.DI using statistical theory has been also defined to evaluate the sensitiveness of SP for condition diagnosis.MPSO algorithm with adaptive inertia weight adjustment and particle mutation is proposed for condition identification.MPSO  ( For a multiresolution of multiplicity  > 1.
Similar to scalar wavelet, Ψ and Φ satisfy the two-scale matrix refinement equations: where  ∈  and  is the set of integers.  and   are lowpass and high-pass matrix filter banks, respectively.
Decomposition GHM multiwavelet constructed by Geronimo, Hardin, and Massopust is one of the most important multiwavelet systems with two pairs of scaling and wavelet functions and has the superior properties of short support, symmetry, orthogonality, and second approximation order [29].Because of the excellent properties, GHM multiwavelets are adopted in this study.The multiscaling functions and multiwavelet functions of GHM multiwavelets are presented in Figure 2. The dilation and wavelet equations for GHM multiwavelet have four coefficients as follows: In view of the matrix filter banks, preprocessing is necessary to translate one stream input signal into two streams.Some preprocessing of preprocessing for multiwavelets has been proposed, such as repeated-row preprocessing, balanced multiwavelet, and prefilter methods [30].Different preprocessing methods will produce different effect on performances of multiwavelets.It is a fundamental problem for each multiwavelet function to choose an appropriate preprocessing method for specific applications.In this study, the preprocessing method of repeated-row preferable for GHM multiwavelet is adopted and given as follows [31]: where   is original signal,   is the signal after preprocessing, and  = √ 2.

Multiwavelet Adaptive Threshold
Denoising.Similar to single wavelet, multiwavelet denoising depends largely on the threshold denoising.The effect of threshold denoising depends on the selection of thresholds.A variety of threshold choosing methods can be mainly divided into two categories: global thresholding and level-dependent thresholding.The former chooses a single value of  to be applied globally to all empirical wavelet coefficients, while the latter chooses different threshold value   for each wavelet level.However, it is difficult to choose appropriate threshold values for different wavelet coefficient.A large threshold value cuts too many coefficients, resulting in the loss of useful information.Conversely, a too small threshold value will leave much noise.In this study, the method of adaptive selecting appropriate threshold values for multiwavelet denoising based on comparison of noise energy in different levels is proposed.
According to the noise levels of wavelet coefficients, the adaptive threshold value is determined by energy ratio.The energy spectrum of multiwavelet coefficient is denoted as follows: where   represents multiwavelet coefficient total energy in the th layer;    represents the multiwavelet coefficient energy of the -dimensional in the th layer;   () is the multiwavelet coefficient in the th layer after -dimensional multiwavelet decomposition; and  is the number of dimensions of multiwavelet coefficient.
Energy ratio of multiwavelet coefficient can be obtained as follows: According to noise level of multiwavelet coefficients, the threshold values of each multiwavelet coefficient can be adaptively obtained as follows: where  is the median absolute value of multiwavelet coefficient;  is signal length.
In conclusion, the processing steps of multiwavelet adaptive threshold denoising are summarized as follows.
(1) Preprocess the original signal to transform it into two streams by the method of repeated-row preferable.
(3) Threshold values are adaptively determined by energy ratio of the wavelet coefficients.
(4) Threshold the wavelet coefficients.(5) Reconstruct the threshold wavelet coefficients and the scale coefficients.
(6) Postprocess the two stream results to get the denoising signal.
In order to test effectiveness of multiwavelet adaptive threshold denoising proposed in this paper, a simulation experiment is designed as follows.
The simulation signal is composed of a periodic impulse component and white Gaussian noise to simulate a bearing fault.The periodic impulse signal with the period of 0.01 s is expressed as where  is damp coefficient;   denotes natural frequency;  0 indicates displacement constant.The shock impulse signal is displayed in Figure 3(a).In this case,  = 0.1,   = 3 kHz,  0 = 5, and sampling frequency and sampling points are 20 kHz and 4096, respectively.The simulation signal is shown in Figure 3(b), the signal has a low signal-noise ratio (SNR), and no useful features can be seen in the dynamic signal in the time domain.The noisy signal is processed using GHM multiwavelet adaptive threshold denoising, GHM multiwavelet neighboring coefficient denoising, and Daubechies 2 (db2) wavelet threshold denoising, respectively.Type of thresholding used is soft thresholding, and decomposition level is four.Denoised signal's performance is evaluated based on mean square error (MSE) and SNR. Figure 4 shows the denoising results using different wavelet denoising techniques.The SNR and MSE of different wavelet denoising techniques are calculated, as shown in Table 1.The results indicate that the method of GHM multiwavelet adaptive threshold denoising has the maximum SNR and the minimum MSE, which means the method proposed in this study can effectively extract Frequency-domain kurtosis: Mean frequency that wave shape cross the mean of timedomain signal: Stabilization factor of wave shape: Sum of the squares of the power spectrum: Square root of the sum of the squares of the power spectrum: where  is the number of spectrum lines,   is frequency, and from 0 Hz to the maximum analysis frequency, (  ) is the power spectrum value at frequency   , and  = The distinction rate (DR) is defined as It is obvious that the larger the value of the DI, the larger the value of the DR will be and, therefore, the better the SP will be.Thus, the DIcan be used as the index of the quality to evaluate the distinguishing sensitivity of the SP.
The number of symptom parameters used for diagnosis and fault types are  and , respectively; the synthetic detection index (SDI) is defined as follows:

MPSO for Condition Diagnosis
4.1.Brief of PSO.PSO algorithm is based on the groups, and according to the environmental fitness, individual in groups will be moved to the good region.The algorithm evaluates the optimal result by using evolutionary fitness function of group, and each particle in the algorithm has a fitness value determined by the fitness function; two properties of position and speed that are used to show the position and moving speed of the current articles in the solving space, by the fitness function value corresponding to particle position coordinate, determine the performance of particles.
In PSO algorithm, each particle adjusts its position according to its own experience and according to the experience of a neighboring particle, making use of the best position encountered by itself and its neighbor.
In the -dimensional search space, the  particle's space position is defined as follows: The velocity of particle  is defined as follows: The best previous position of particle  is defined as follows: The best position among all particles experienced is defined as follows: The particle updates the position and velocity according to the following equations: where  1 and  2 are the random numbers within (0, 1) and  1 and  2 are the acceleration which constants the control of how far a particle moves in a single generation.The inertia weight  controls the previous velocity of particle, and it is defined as follows: where rand is random generated number between 0 and 1.
Although PSO algorithm is easy to realize, the method is easy to trap into local optimum.Shi and Eberhart proposed a linearly decreasing weight particle swarm optimization (WPSO) of which a linearly decreasing inertia factor was introduced into the velocity of the updated equation from the original PSO [35,36].The performance of WPSO is significantly improved over the original PSO because WPSO balances out the global and local search abilities of the swarm effectively.The equation for the linearly decreased weight is defined as follows: where  max is 1,  min is 0.1, and iteration max is the maximum number of the allowed iterations.
The velocity of the updated equation for WPSO is defined as follows:

MPSO.
Although WPSO algorithm improved conventional PSO to a certain extent, it cannot adapt to all of complex practical problems.The main reasons can be explained as follows.(1) The inertia weight of conventional WPSO algorithm is monotone decreasing, and adjustment ability of WPSO algorithm is limited.If particles cannot find optimal point in the initial stage of the algorithm, WPSO algorithm is easy to trap into local optimum with the decrease of the inertia weight.( 2) With increasing iterations, particle diversity of WPSO algorithm decreases; it causes deterioration of global search ability; WPSO algorithm is also easy to trap into local optimum and premature convergence.
To improve global search ability and adjustment ability of conventional PSO algorithm and prevent local optimum and premature convergence problems, MPSO algorithm with adaptive inertia weight adjustment and particle mutation is proposed in this paper.
where () is optimum fitness value of the th iteration; ( + 5) is optimum fitness value of the ( + 5)th iteration;  indicates change rate of fitness value in five iterations.
According to the variation of , the inertia weight  adaptively adjusts as follows: where  is a random number with a uniform probability within 0∼1;  1 and  2 are parameters;  1 >  2 ; the choice of  1 and  2 is determined experimentally; here  1 = 0.5 and  2 = 0.2.When  > 0.05, the algorithm is in the exploration stage, and a large  is beneficial to the algorithm's convergence.
When  ≤ 0.05, the algorithm is in the development stage, and a small  is beneficial to searching optimum point.

Particle Mutation.
To increase particle diversity of PSO algorithm, the method of particle mutation is proposed.In the operation process of PSO algorithm, if the best position among all particles best does not change in a long time, some particles are mutated according to a certain probability.The execution process of the mutation for PSO is as follows.
(1) All particles are arranged in ascending order according to the values of the fitness function.(2) The  ( > 1) particles with smaller fitness functions are selected.
(5)   is compared with   , if   >   , and then the particle's space position is updated by using ( 29). ( 6) Steps ( 3)-( 5) are looped until the space position of  particles are updated: where  is random data that obeys Gaussian(0, 1) distribution.

Fitness Function of MPSO for Condition Diagnosis.
Assume that  is the sample set of vibration signals measured in  different states; the length of  is ,  = { 1 ,  2 , . . .,   }.
Every sample signal has  identified symptoms (in this paper, the symptoms are  1 - 6 ).Then, the clustering analysis is to divide  sample data into  states, such that the fitness function  shown in (30) is minimized:

Diagnosis and Application
In this section, the application of condition diagnosis for a rolling bearing is shown to verify that the method proposed in this paper is effective.

Diagnosis by the Proposed Method.
The main procedure for fault diagnosis using GHM multiwavelet adaptive threshold denoising and MPSO algorithm is shown in Figure 11 and explained as follows.
(1) Vibration signals are measured in each known state.
(2) Weak fault feature is extracted by using GHM multiwavelet adaptive threshold denoising.
(4) The highly sensitive SPs are selected for condition diagnosis by DI.
(5) MPSO algorithm is trained with SPs selected by DI, and the optimal clustering centers are obtained.In this study, the good SPs which have high sensitivity for distinguishing each fault state of the bearing are selected by the method of DI.As an example, Table 2 lists parts of DIs of SPs.The maximum value (50.49) of SDI is obtained in the case of the combination of  1 ,  5 , and  6 , and when  1 ,  5 , and  6 , are used for distinguishing each state separately, the DIs are larger than 1.41, and all of the DRs are larger than 92.1%.Therefore,  1 ,  5 , and  6 have high sensitivity for distinguishing each fault state of the bearing.
In this study, MPSO automatically obtains the optimal clustering centers according to the classification of the sample data information.The purpose of training MPSO is the acquisition of optimum clustering centers.The SPs selected by DI were input into MPSO.MPSO converged to the optimum clustering centers.In the training process of MPSO, at first, the sample data are classified into the normal, the outer-race defect, the inner-race defect, and the roller element defect randomly.The fitness values and the clustering centers are calculated by (30) and (31).With increasing iterations, the speed and position (classification of the sample data) of the particle are updated incessantly, and according to the classification of the sample data information, the clustering centers are also updated.Finally, the optimal clustering centers with a minimum fitness value are calculated.
To explain the effectiveness of MPSO algorithm, a comparison is made among MPSO, WPSO, and PSO algorithms.The optimal clustering centers of each state are obtained by MPSO, WPSO, and PSO algorithms, respectively.Particle number and iteration number of the three algorithms are 50 and 1000, respectively.The clustering centers obtained by each method are shown in Tables 3, 4, and 5. Figure 12 shows the fitness curve of PSO, WPSO, and MPSO algorithms.It is obvious that MPSO algorithm has the minimum fitness value; namely, the optimal clustering centers obtained by MPSO algorithm are the most accurate.MPSO algorithm has stronger capacity of searching optimal solution.
After training MPSO algorithm, to verify the diagnostic capability of the proposed method in this paper, the test data measured in each known state that had not been used to train MPSO algorithm were used.When inputting the test data into the trained MPSO algorithm, MPSO classified the test data according to the information of the optimum clustering centers shown in Table 3 and correctly and quickly output identification results.As an example, some diagnosis results are listed in Table 6.We also identified condition of the rolling bearing using PSO and WPSO algorithms, respectively, and   7 and 8.The comparison of diagnostic capability of each method is shown in Figure 13.Viewing the overall diagnostic results, diagnostic accuracy of each state using MPSO algorithm is 100%, 88.3%, 86.7%, and 81.7%, respectively; they are the largest in three methods.The method proposed in this study provides a more accurate estimate in the case of the rolling bearing faults diagnosis.These results verified the efficiency of the intelligent diagnosis method using multiwavelet adaptive threshold denoising and MPSO proposed in this paper.

Conclusions
In order to diagnose faults of rotation machinery at an early stage, this paper proposed a novel intelligent condition diagnosis method using multiwavelet adaptive threshold denoising and MPSO to detect faults and distinguish fault types at an early stage.The main conclusions are summarized as follows.
(1) The method of multiwavelet adaptive threshold denoising was presented for extracting weak fault features under background noise.It could adaptively select appropriate threshold for multiwavelet with energy ratio of multiwavelet coefficient.The simulation experiment verified that the method of multiwavelet adaptive threshold denoising can effectively extract fault features and eliminate much noise from the noisy signal.
(2) The six SPs in the frequency domain were defined for reflecting the features of vibration signals measured in each state.DI using statistical theory had been also defined to evaluate the applicability of the SPs for the condition diagnosis measured in each state.DI could (3) MPSO algorithm with adaptive inertia weight adjustment and particle mutation was proposed for condition identification.MPSO algorithm was used to classify the SPs calculated from the signals in each machine state for condition diagnosis, as well as obtaining their optimal clustering centers.According to these optimal clustering centers' information, the conditions of rotation machinery could be accurately identified.MPSO algorithm effectively solved local

Figure 3 :Figure 4 :
Figure 3: The simulation signal in the time domain: (a) the shock impulse signal; (b) the noisy signal.

Figure 5 Figure 5 :
Figure 5: Experimental system for bearing fault diagnosis.

( 6 )
Condition of the bearing can be diagnosed by the trained MPSO algorithm and SPs.

Figure 8 :Figure 9 :
Figure 8: The vibration signal in outer-race defect state: (a) original vibration signal, (b) after multiwavelet adaptive threshold denoising, and (c) Fourier spectrum of the denoised signal.

Figure 10 :
Figure 10: The vibration signal in roller defect state: (a) original vibration signal, (b) after multiwavelet adaptive threshold denoising, and (c) Fourier spectrum of the denoised signal.

Figure 11 :
Figure 11: Flowchart for the condition diagnosis by the method proposed in this study.

Table 1 :
SNR and MSE of different wavelet denoising techniques.
1 ∼ . is mean value of the analysis frequency, and  = (∑  =1   ⋅ (  ))/ ∑  =1 (  );  is standard deviation, and  = √ (∑  =1 (  − ) 2 ⋅ (  ))/.3.2.Detection Index.For automatic diagnosis, SPs are needed that can sensitively distinguish the fault types.In order to evaluate the sensitivity of a SP for distinguishing two states, such as a normal or an abnormal state, DI is defined as follows.Supposing that  1 and  2 are the SP values calculated from the signals measured in state 1 and state 2, respectively, their average value and standard deviation are  and .The DI is calculated by DI

Table 2 :
DIs of each SP.

Table 4 :
Clustering centers obtained by SPO.
some identification results are shown in Tables