A Rolling Element Bearing Fault Diagnosis Approach Based on Multifractal Theory and Gray Relation Theory

Bearing failure is one of the dominant causes of failure and breakdowns in rotating machinery, leading to huge economic loss. Aiming at the nonstationary and nonlinear characteristics of bearing vibration signals as well as the complexity of condition-indicating information distribution in the signals, a novel rolling element bearing fault diagnosis method based on multifractal theory and gray relation theory was proposed in the paper. Firstly, a generalized multifractal dimension algorithm was developed to extract the characteristic vectors of fault features from the bearing vibration signals, which can offer more meaningful and distinguishing information reflecting different bearing health status in comparison with conventional single fractal dimension. After feature extraction by multifractal dimensions, an adaptive gray relation algorithm was applied to implement an automated bearing fault pattern recognition. The experimental results show that the proposed method can identify various bearing fault types as well as severities effectively and accurately.


Introduction
As a main component in almost all types of rotating machine rolling element bearings have been used widely. Bearing failure is one of the main causes of failure and breakdowns in rotating machinery, leading to significant economic loss [1][2][3]. So as to keep machinery running with high reliability and availability, it is necessary to develop a reliable and effective bearing fault diagnosis method. Among various bearing fault diagnosis methods, vibration-based analysis methods have arrested extensive attention in the last few decades [4]. Bearing vibration signals contain rich information about machine health status, and thus it is possible to obtain foremost characteristics from vibration signals by signal processing techniques [5]. Nowadays a great number of signal processing methods have been used for rolling bearing fault detection and diagnosis. However, due to the nonlinear factors, e.g., friction, stiffness and clearance, bearing vibration signals in particular under faulty conditions, will behave nonstationary and nonlinear characteristics [6]. In addition, the measured vibration signals contain not only the working condition information concerning the rolling bearing itself, but abundant a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 information of other moving parts and structures of machinery, which belongs to the category of background noise compared with the former [7]. As background noise is often comparatively large, the slight bearing fault information may easily be submerged in the background noise and become hard to extract and detect. Therefore, the conventional time and frequency domain methods which aim mainly at linear vibration signals, and even advanced signal processing methods (e.g., wavelet transform) may not make an accurate assessment of the rolling bearing health status [8].
With the development of nonlinear dynamics theory, a great deal of nonlinear analysis techniques have been proposed aiming at recognizing and predicting the complex nonlinear dynamic behavior of bearings [4,9]. One of the most common manners is to refine and extract the fault features from vibration signals by combination of a few of advanced signal processing methods (e.g., wavelet package transform (WPT) [8,10], hilbert transform(HT) [10,11], empirical mode decomposition (EMD) [11] and higher order spectra (HOS) [12]), to recognize the fault frequency used to compare with the theoretical value with involvement of expert's empirical judgement. With the advent of artificial intelligence, the process of rolling bearing fault diagnosis is more and more treated as a procedure of fault pattern recognition, and its effectiveness and reliability significantly depend on the selection of dominant characteristic vector of the fault features [13]. Recently, Some entropy based methods, e.g., approximate entropy (ApEn) [14,15], sample entropy (SampEn) [16], fuzzy entropy (FuzzyEn) [16,17], hierarchical entropy (HE) [13,18] and hierarchical fuzzy entropy (HFE) [13], have been proposed to to extract characteristic vector of the fault features from bearing vibration signals. In the paper, we introduce a multifractal theory based method, i.e., a generalized multifractal dimension algorithm, to extract dominant characteristic vector of the fault features from bearing vibration signals. Fractal theory is one of the most important branches for the contemporary nonlinear science, and is in particular suitable for processing all kinds of complex nonlinear and nonstationary phenomenon [19] and thus may also suitable for fault feature extraction from bearing vibration signals.
Usually, after feature extraction, a pattern recognition technique is needed to achieve the rolling element bearing fault diagnosis automatically [13]. A series of pattern recognition methods have been introduced for mechanical fault diagnosis nowadays, among which the most commonly used ones are artificial neural networks (ANNs) [20,21,22] and support vector machines (SVMs) [23]. The training of ANNs always needs a large number of samples which may be difficult or even impossible to achieve in practical applications, especially for the fault ones. While SVMs based on statistical learning theory, which is of specialties for a smaller number of samples, have better generalization than ANNs and ensure that the local and global optimal solution are exactly the same [24]. However, the accuracy of a support vector machine (SVM) classifier is highly decided by the selection of optimal parameters for SVMs [24,25]. To ensure the diagnostic accuracy, an optimization algorithm [24,25] or / and complex multiclass concept [13,26] have to be used subsidiarily to improve the effective of SVMs. Here, in order to solve the issue of generality versus accuracy, an adaptive gray relation algorithm was developed to achieve accurate patter recognition based on a small number of samples.
In summary, a novel approach for rolling element bearing fault diagnosis was proposed based on multifractal theory and gray relation theory. First, the fault features from bearing vibration signals, which offer more useful and distinguishing information imaging different bearing health status were extracted by generalized multifractal dimension algorithm. And then different fault types of rolling element bearing as well as severities were recognized by an adaptive gray relation algorithm. The rest of this paper is organized as follows. A review of conventional single fractal dimension and its improved version is introduced separately in section 2.1 and 2.2. The gray relation theory for fault pattern recognition is given briefly in section 3, followed by section 4 which presents the proposed bearing fault diagnosis method. Next, experimental validation of the proposed method is given in section 5.1 and 5.2. Finally, conclusions are illustrated in section 6.

Fractal theory
Fractal theory is one of the most important branches for the contemporary nonlinear science, and is a novel scientific method especially suitable for processing various complex, nonlinear, irregular and nonstationary phenomena [27], and thereafter may be also suitable for fault feature extraction from bearing vibration signals.

Traditional single fractal dimension
Fractal theory has been used to radar signal analysis [27], in which fractal dimensions such as box-counting dimension, Higuchi dimension, Petrosian dimension, Katz dimension, Sevcik dimension and so forth [28], are applied to depict the complexity of radar emitter signals. And among them, box-counting dimension is one of the most common fractal dimensions due to simplicity of its algorithm, and its algorithm is described as follows.
Given that A is a nonempty bounded subset of Euclidean space R n to be calculated, and N (A,ε) is the least number of boxes covering A with the side length of ε, and therefore the boxcounting dimension can be defined as: For the actual bearing vibration discrete signals, due to existence of the sampling frequency, the sampling interval σ is the highest resolution for the signals. Therefore, it is no meaning to calculate the box-counting dimension when ε ! 0 and the minimum length of the box is often taken as ε = σ.
Considering a sequence of vibration discrete signals x(i) as a closed set of Euclidean space R n , the computational procedure of box-counting dimension is described as follows: Adopt approximate method to make the minimum side length of the boxes covering the sequence of vibration discrete signals x(i) equal to the sampling interval σ, and then calculate the least number of boxes N kσ with side length of kσ covering the sequence of vibration discrete signals x(i), thus: pðkεÞ ¼ where i = 1,2,Á Á Á,N 0 /k, k = 1,2,. . .,K. N 0 is the number of sampling points, K < N 0 . p(kε) is the longitudinal coordinate scale of the sequence of vibration discrete signals x(i). Therefore N kε can be expressed as: Select a good linearity of the fitting curve log kε * log N kε as the scale-free area, thus: where k 1 k k 2 , k 1 and k 2 are the beginning and end of the scale-free area respectively.
Usually, the least square method can be used to calculate the slope of the current fitting curve which is the fractal box-counting dimension D of the sequence of vibration discrete signals x(i): However, single fractal dimension is often not enough to describe a complicated fractal object, and thus multifractal theory is introduced to analyze bearing vibration signals. Multifractal dimensions are the extension of conventional single fractal dimension, which can be applied to describe growth features of a fractal object at different levels, to compensate the lack of conventional single fractal dimension.

Generalized multifractal dimensions
The traditional single fractal dimension has been introduced in section 2.1, and has been widespreadly used in strictly self-similar signals such as the signals in the biological medicine, image analysis and electromagnetic fault diagnosis. Nowadays, single fractal dimension has also been used to quantitatively extract fault features from bearing vibration signals, which has attracted extensive attention from investigators. However, the common bearing vibration signals do not satisfy the self-similar structure of fractal theory to some degree. Therefore, when using the fractal box-counting dimension to calculate box-counting dimension of the vibration signals, the curve fitting often do not have good linear structure. Moreover, single fractal dimension only images the whole characteristics of vibration signals and lacks the depiction of local singularity of the vibration signals. In comparison, Multifractal dimensions are the extension of conventional single fractal dimension, which can be applied to describe growth features of the fractal object at various levels, to compensate the lack of conventional single fractal dimension. And the detail is described as follows.
Divide the research object (i.e., bearing vibration signals) into N small areas, and consider the linear dimension of the i th area is ε i , and thus the probability density function P i for the i th area with different scale index α i can be described as: Define : where X q (ε) is the weight sum of probability of each area; q is the power of the probability density function P i for the i th area.
where D q is the generalized multifractal dimensions. From Eq (10), it can be seen that when q>>1, the areas with bigger probability play a major role in X q (ε), and X q (ε) and D q mirror the characteristics of areas with bigger probability (i.e., areas with denseness). On the contrary, when q<<1, X q (ε) and D q mirror the characteristics of areas with smaller probability (i.e., areas with sparseness). Therefore, the characteristics of different areas with different probability are reflected by different q, and a sequence of signals are divided into many areas, which are full of different singularity after weight sum of probability of each area. And thus the subtle structure of a sequence of signals can be described at different levels by the generalized multifractal dimensions.

Gray relation theory
After fault feature extraction, a fault pattern recognition technique is used to achieve automatically the rolling element bearing fault diagnosis. The study on gray relation theory is the foundation of gray system theory, which is based on the basic theory of space mathematics to calculate relation coefficient and relation degree between reference characteristic vectors and comparative characteristic vectors. Many investigations [29][30][31][32] have demonstrated that gray relation theory is full of capability to be used in machinery fault recognition with four reasons: 1. it has the ability to assist the selection of characteristic parameters for classification; 2. it has good tolerance to measurement noise; 3. it can solve the learning problem with a small number of samples; 4. its algorithm is simple and can alleviate the issue of generality versus accuracy.
Given the dominant characteristic vectors (i.e., the multifractal dimensions) of the fault features extracted from bearing vibration signals to be recognized is as follows.
where B i (i = 1,2,. . .) is fault pattern to be identified; D k,i (i = 1,2,. . .) is each characteristic parameter; K is the total number of characteristic parameters selected as characteristic vector. Consider the knowledge base from a small number of samples between fault signatures (i.e., the characteristic vectors) and fault patterns (i.e., the fault types of rolling element bearings as well as various levels of severity) is as follows.
where ρ is distinguishing coefficient, and its value is usually set as 0.5; ξ(D k,i ,c j (k)) is the gray relation coefficient for k th characteristic parameter of B i and C j ; ξ(B i ,C j ) is the gray relation degree of B i and C j . And then B i is classified to the fault pattern in which the maximum ξ(B i ,C j )(j = 1,2,. . .,) is obtained.
So as to enhance its tolerance to measurement noise and the ability to assist the selection of characteristic values for fault classification, the information theory was introduced into the calculation of the relation degree and that is so-called adaptive gray relation algorithm (GRA) [33].
First, process the distance of characteristic parameter |Δd ij (k)| = |D k,i −c j (k)| and then calculate the probability as follows.
where M is the total number of the known fault patterns in the knowledge base. Define the entropy value as follows.
And the maximum entropy value is as follows.
Then, the relative entropy value is obtained as follows.
Referred to the concept of surplus degree in information theory, the definition of surplus degree for k th characteristic parameter is as follows.
The essential meaning of surplus degree is to remove the difference between the entropy value of the k th characteristic parameter and the optimal entropy value of characteristic parameter. The bigger H ij (k) is, the more important the k th characteristic parameter is, and then the H ij (k) should be granted greater weight.
Finally, calculate the weight coefficient a ij (k) for k th characteristic parameter as follows.
Then obtain relation degree by weight coefficient multiplying with the corresponding relation coefficient as follows.
And finally B i is classified to the fault pattern in which the maximum ξ(B i ,C j )(j = 1,2,. . .,M) is obtained, i.e., the the fault types of rolling element bearings as well as various levels of severity are recognized.

Proposed Approach
The extraction of subtle features from bearing vibration signals is the core of a rolling element bearing fault diagnosis, and its quality will directly produce essential effect on the success rate of the subsequent fault pattern recognition. With the integration of the advantages of multifractal dimension algorithm and adaptive gray relation algorithm, a novel fault diagnosis method for rolling element bearings is put forward as follows.
1. The bearing vibration signals under different health status are acquired as samples which are divided into two subsets, the samples for knowledge base and the samples for testing.
2. The dominant characteristic vectors of fault features from the bearing vibration signals, which can provide more useful information concerning bearing health status, were extracted by multifractal dimension algorithm.
3. Establish knowledge base according to the relationship between fault signatures (i.e., the characteristic vectors) and fault patterns (i.e., the fault types of rolling element bearings as well as different levels of severity) from the base samples for adaptive GRA model. 4. The feature vectors of the testing samples are input into the adaptive GRA model, and various bearing health status can be recognized by the output of the adaptive GRA model.
And the diagnostic procedure is shown in Fig 1. The detailed procedure of the proposed method is as follows: 1) Discretize the unknown received bearing vibration signals: Suppose the received bearing vibration signals is s, and the sequence of preprocessed discrete bearing vibration signals is {s(i)}, i = 1,2,Á Á Á,N 0 , where N 0 is the length of the sequence.
2) Recombine the sequence of discrete bearing vibration signals {s(i)}, i = 1,2,Á Á Á,N 0 : where n is the time of recombining the sequence of discrete bearing vibration signals.
3) Extract multifractal dimensions from the total sequence of recombined discrete bearing vibration signals S(j), j = 1,2,Á Á Á,n: Divide the research object (i.e., the total sequences of recombined discrete bearing vibration signals S(j), j = 1,2,Á Á Á,n) into N number of small areas respectively, and consider the linear dimension of the i th area is ε i , and thus the probability density function P i for the i th area with different scale index α i can be described as: where α i is singular index and is the fractal dimension for i th area. As the total sequences of recombined discrete bearing vibration signals S(j), j = 1,2,Á Á Á,n, has been divided into N number of small areas, thus a series of α i (i = 1,2,Á Á Á,N) can be obtained and made up of a variable f (α) which is the multifractal spectrum for each sequence of recombined discrete bearing vibration signals.
Define : X q ðεÞ ¼ where X q (ε) is the weight sum of probability of each area; q is the power of the probability density function P i for the i th area.
where D q is a generalized multifractal dimension. Obtain the sum S j of the j th sequence of recombined discrete bearing vibration signals S(j): Obtain the sum S of the whole discrete bearing vibration signals {s(i)}, i = 1,2,. . .,N 0 : And then calculate probability density function P j for the j th area: The generalized multifractal dimension D q can be obtained by inputting Eq (31) to Eq (28). 4) Establish knowledge base for adaptive GRA model and various unknown bearing health status can be recognized by the output of the adaptive GRA model: Set q as from −q 0 to q 0 , and then there are 2q 0 + 1 sets of multifractal dimensions, in which there are n = log 2 N 0 number of characteristic points for each q. Thus the total number of characteristic parameters for a sequence of preprocessed discrete bearing vibration signals {s(i)}, i = 1,2,. . .,N 0 , is M = (2q 0 + 1)Án = (2q 0 + 1)Álog 2 N 0 , to constitute one set of characteristic vector. And finally the unknown bearing health status can be classified to the fault pattern in which the maximum relation degree is obtained, i.e., the the fault types of rolling element bearing as well as various levels of severity is recognized by Eq (15) to Eq (21).

Experimental rig
All the rolling element bearing vibration signals used for analysis were downloaded from Case Western Reserve University Bearing Data Center [34]. The deep groove balling bearing 6205-2RS JEM SKF was used in the experiment. The whole experimental rig consists of a two horsepower three-phase induction motor, a torque transducer and a dynamometer, as shown in Fig 2. The desired torque load levels could be achieved by controlling the dynamometer. The horsepower and speed data are collected by the sensor. The motor shaft at the drive end is supported by the test bearing, and an accelerometer with a bandwidth up to 5000Hz was installed on the motor housing at the drive end of the motor, and then bearing vibration data under different working conditions were collected using a recorder with a sampling frequency of 12 kHz. Here single point faults were introduced into the test bearing by using electro-discharge

Application and analysis
To assess the effectiveness of the proposed method, experimental analyses of rolling element bearing fault diagnosis were carried out. The bearing vibration signals used for the analyses were obtained from the tests which were conducted at the load of 0 horsepower and at the motor speed of 1797 r/min. Normal and three fault types of bearing vibration data as well as those with various severities for each fault type are analyzed, and detailed description of the related vibration data was shown in Table 1. Considering various fault categories and severities, the rolling element bearing fault diagnosis turns out to be an 11-class fault pattern recognition problem. The data set was made up of totally 550 data samples, in which each data sample was cut into a 2048-point time series from the original vibration data and there is no overlap between any two of them. Among these 550 data samples, 110 samples were selected at random as the samples for knowledge base, and the rest 440 were automatically treated as testing data. Here the testing data are four times larger than the amount of the training data, as the bearing vibration data under faulty conditions is hard to obtain and is of small sample size as usual in the practical applications.
The fault features extracted from bearing normal condition and various fault conditions with fault size 7mils (seen Fig 3) by traditional single fractal dimension (i.e., box-counting dimension) was shown in Table 2 and by generalized multifractal dimensions were shown in Fig 4. And the fault features extracted from bearing inner race fault condition with different levels of severity (seen Fig 5) by single fractal dimension (i.e., box-counting dimension) was shown in Table 3  After establishing knowledge base according to the relationship between fault signatures (i.e., the characteristic vectors) and fault patterns (i.e., the fault types of rolling element bearings as well as different levels of severity) from the base samples for adaptive GRA model, the characteristic vectors of the testing samples were input into the adaptive GRA model, and various bearing health status were recognized by the output of the adaptive GRA model shown in Table 4.
From Table 4, the identification results show that the recognition effect with traditional single fractal dimension (i.e., box-counting dimension) is poor and misleading due to the limited fault features extracted, while the recognition effect with generalized multifractal dimensions seems much better, in which the recognition success rate for detecting bearing faulty conditions is 100% and the total success rate is 96.59% for identifying different bearing fault types as well as severities.

Further discussion
To further analyze the effectiveness of the proposed method for dealing with the learning problem with an extremely small number of samples, another experimental test in which only one random base sample for these 11 classifications used for establishing knowledge base for the adaptive GRA model was carried out. The feature vectors of the testing samples were also input into the adaptive GRA model, and various bearing health status were recognized by the output of the adaptive GRA model shown in Table 5.
From Table 5, it is encouraging to see that the proposed method shows good robust recognition effect over this extreme experimental validation and the total success rate can still reach more than 81%, while the recognition success rate for detecting bearing faulty conditions is still 100%.
In the future research, for improving the diagnostic accuracy of the proposed approach, other advanced signal processing methods (e.g., wavelet package transform (WPT), hilbert  transform (HT), empirical mode decomposition (EMD) or higher order spectra (HOS)), should be explored to be integrated into the generalized multifractal dimensions to more effectively extract dominant characteristic vectors.

Conclusion
Rolling element bearings as an important component in almost all types of rotating machines have been widespreadly used and its failure is one of the foremost causes of failure and breakdowns in rotating machinery, resulting in significant economic loss. In the paper, a novel approach for rolling element bearing fault diagnosis was proposed based on generalized multifractal dimension algorithm and adaptive gray relation algorithm. First, the fault features from the bearing vibration signals, which can provide more useful information reflecting bearing health status were extracted by the generalized multifractal dimension algorithm. And then, the fault types of rolling element bearings as well as different levels of severity are recognized by the outputs of the adaptive GRA algorithm. The experimental results demonstrate that the proposed method can effectively and accurately identify different bearing fault types as well as severities. And some other meaningful conclusions can be obtained as follows: 1. Traditional single fractal dimension is not enough to describe the rolling element bearing health status, while the multifractal dimensions are suitable to be applied to extract fault features from the rolling element bearing vibration signals, which offer more meaningful and distinguishing information reflecting different bearing health status.   2. Gray relation theory is full of capability to be used in the rolling element bearing fault pattern recognition, and the recognition success rate for detecting bearing faulty condtions is 100% and the total success rate can be more than 96% for identifying different bearing fault types as well as severities by the proposed method.
3. When the number of samples for knowledge base becomes extremely limited, the total success rate may be reduced significantly for identifying different bearing fault types as well as severities, while the recognition success rate for detecting bearing faulty conditions is still 100%.