Application of Multi-Scale PCA and Energy Spectrum to Bearing Fault Analysis and Detection in Rotating Machinery Текст научной статьи по специальности «Медицинские технологии»

Application of Multi-Scale PCA and Energy Spectrum to Bearing Fault Analysis and Detection

in Rotating Machinery

K. Baiche, M. Zelmat, A. Lachouri

Abstract — This paper introduces various works on fault detection of rotating machines caused by bearings damage. The novelty of this work is the application of new scheme based on the multi-scale principal component analysis MSPCA, and the energy spectrum of details coefficients of the discrete wavelet transform (DWT). The DWT coefficient details are calculate in first stage in order to be used as inputs of the MSPCA scheme and in the second stage they are used to evaluate the spectrum energy.

Keywords: Fault analysis, rotating machinery, MSPCA, wavelet transform, spectrum energy, bearing diagnosis.

I. Introduction

Rolling element bearings are essential elements in most rotating machines. Bearing failures are prone to cause both personal damage and economic loss if not be detected well in advance [1-5]. To avoid such catastrophic failures, it is necessary to develop and implement efficient diagnosis monitoring systems that are independent of operating conditions. Therefore, a significant amount of research efforts have focused on the predictive maintenance of machines. Machine Vibration Signature Analysis (MVSA) provides an important way to assess the health of a machine.

In traditional MVSA, the Fourier transform is used to determine the vibration spectrum. Typically, frequencies of bearing defects are identified and compared with initial measurements to detect any deterioration in bearing health.

In recent years [5], different technologies have been used in order to process signals produced by dynamical systems. Most of the authors classify the analysis of vibration signature in three approaches. The first one is the time domain, based on statistical parameters such as mean, root mean-square, variance, kurtosis, etc. The second approach proceeds in the frequency domain, where the Fourier transform and its numerous variants have been intensively used. The shortcoming of this approach is that Fourier analysis is theoretically limited to stationary signals, whilst bearing vibrations are typically non-stationary by nature

[6]. The third approach is based on the time-frequency analysis such as the short-time Fourier transform (STFT) or the (discrete) wavelet transform (DWT). It has been successfully applied as a fault feature extractor due to its

Manuscript received November 18, 2012.

K. Baiche is with Applied Automatic Laboratory, LAA, University of Boumerdes, Algeria (e-mail: [email protected]).

M. Zelmat is with LAA, University of Boumerdes, Algeria (e-mail: [email protected]).

A. Lachouri is with LAA, University of Skikda, Algeria (e-mail: [email protected]).

good energy concentration properties. Besides, in the diagnosis field, methods based on the information redundancy concept have been developed. Generally, these methods rest on a consistency test between an observed behavior of the process provided by sensors and a mathematical representation. The comparison between these behaviors produces residual quantities that serve to discriminate normal plant operation from an abnormal situation.

In this work, experimental vibration signals for both normal and faulty bearings are acquired. Then a feature vector is generated to be applied as the inputs to a MSPCA according the processing flowchart in Fig.1. The square prediction error (SPE) calculated from the MSPCA method and the spectrum Energy (SE) indicate the difference between the healthy and the defective bearing.

Coef Di

Fig. 1. Processing flowchart

This manuscript develops two methods based on the spectrum energy and the MSCPA. Principal components analysis (PCA) is a projection based on statistical methods used for dimensionality reduction. It is a descriptive technique which permits the study of the relations existing between variables without taking account of any a priori structure [7].

The discrete wavelet transform (DWT) provides an efficient method for generating feature vectors. The DWT coefficients can be used to generate statistical parameters from each resolution level of the transform. These coefficients capture trends within the sensors at the corresponding scale [8]. They serve as efficient inputs to be fed to a decision system, such a Neural Network.

The remainder of this paper is organized as follows. The next section describes the characteristic vibration frequencies of bearings. Section III presents the experiment setup and section IV describes the PCA, wavelets and MSPCA methods. The Section V describes the spectrum Energy and the Section VI presents the simulation, and experimental results, respectively. The conclusion is given in section VII.

R&I, 2012, №4

31

II. Characteristic vibration frequencies

The characteristic vibration frequencies due to bearing defects can be calculated from the rotor speed and bearing geometry. The typical rolling element bearing geometry is displayed in Fig.2. Characteristic vibration frequency fcan be calculated using (1)-(4) [1], [2].

The outer race defect frequency, fOD -the ball passing frequency on the outer race- is given by,

fOD = - frm ^ - YD c°s^| (1)

where q is the contact angle, PD is the pitch diameter, BD is the ball diameter, n is the number of balls, and f rm is the rotational speed.

Pitch diameter PD

Fig. 2. Rolling element bearing geometry

discharge machining with fault diameters of 7 mils (1 mil=0.001 inches) and depth of 110 mils in inner race. SKF bearings were used. Drive end (6205-2RS JEM SKF, deep groove ball bearing) and fan end (6203-2RS JEM SKF, deep groove ball bearing) bearing specifications, including bearing geometry and defect frequencies are listed in Tables I-IV. Vibration data was collected using accelerometers, which were attached to the housing with magnetic bases.

Fig. 3. Test stand

TABLE I

Drive end bearing geometry (size in inches)

Inside Outside Thickness Ball Pitch Diameter Diameter Diameter Diameter

The inner race defect frequency, fD-the ball passing frequency on the inner race- is given by,

fD

-f

2 r

, BD

i +--cosq

PD

(2)

The ball defect frequency, fBD - the ball spin frequency -is given by

fBD

PD

2BD

f |1

J rm

(3)

The train defect frequency, fD-caused by an irregularity in the train- is given by,

fTD = 1 f rm ^ - jDCOSP) . (4)

The characteristic current frequencies, fCF-due to the bearing characteristic vibration frequencies, fv- are calculated by

fCF = | fe ± mfv\ (5)

where m=1, 2, 3, and feis the power line frequency. This latter equation represents an amplitude modulation of the stator current by the bearing vibrations. Equations (1)-(5) can be used to calculate the current spectral components due to faulty bearings.

III. Description of the experiment

As shown in Fig.3, the test stand consists of a 2 HP motor (left), a torque transducer/encoder (center), a dynamometer (right), and control electronics (not shown). The test bearings support the motor shaft. Single point faults were introduced to the test bearings using electro-

0.9843 2.0472 0.5906 0.3126 1.537

TABLE II

Defect frequencies (multiple of running speed in Hz) for Drive end bearing

Inner Ring Outer Ring Cage Train Rolling Element

5.4152 3.5848 0.39828 4.7135

TABLE III

______FAN end bearing geometry (size in inches)

Inside Outside Thickness Ball Pitch Diameter Diameter Diameter Diameter

0.6693 1.5748 0.4724 0.2656 1.122

TABLE IV

Defect frequencies (multiple of running speed in Hz)

_______________FOR FAN END BEARING______________

Inner Ring Outer Ring Cage Train Rolling Element

4.9469 3.0530 0.3817 3.9874

Accelerometers were placed at the 12 o’clock position at both the drive end and fan end of the motor housing. During some experiments, an accelerometer was attached to the motor supporting base plate as well. Vibration signals were collected using a 16 channel DAT recorder. Digital data was collected at 12,000 samples per second, and data was also collected at 48,000 samples per second for drive end bearing faults. Speed and horsepower data were

32

R&I, 2012, №4

collected using the torque transducer/encoder and were recorded by hand. Outer raceway faults are stationary faults; therefore placement of the fault relative to the load zone of the bearing has a direct impact on the vibration response of the motor/bearing system. In order to quantify this effect, experiments were conducted for both fan and drive end bearings with outer raceway faults located at 3 o’clock (directly in the load zone), at 6 o’clock (orthogonal to the load zone), and at 12 o’clock. The rotation speed (frequency) of the shaft is approximately 1750 rpm (29,2 Hz).

IV. Principal component analysis and wavelets

Generally, the fault detection based on analytic model corresponds to the generation of a residual between normal plant operation and available measurements. From the residual analysis, we can make a decision to indicate whether the failure is present or not.

In this section, the general principles of using PCA for fault detection are presented. it is followed by a brief introduction to wavelets. This section prepares the basic concept of the MSPCA which is explained in the next section.

4.1 Principal component analysis

Principal component analysis (PCA) has been introduced by Pearson and developed by Hotelling. It is used in various domains (meteorology, economy, biology...). Principal component analysis is a descriptive technique; it permits the study of the relations between variables without taking account of any a priori structure. It may also be viewed as a projection tool traditionally used for dimensionality reduction. Consider a data matrix

x(k) = [,x2,l ,xm ]T еЭТm consisting of n sample rows and m variable columns that are normalized to zero mean and unit variance [5].

Once a PCA model is built, a new data sample x, is to be tested for fault detection. It is first scaled and then decomposed as follows,

x = X + x (6)

where

x = PPTx e SP (7)

Is the projection on the principal component subspace (PCS), SP, and,

X = (I - ppt )x e Sr (8)

Is the projection on the residual subspace (Rs), Sr .

For fault detection in the new sample x, a deviation in x from the normal correlation would change the projection into the subspaces, either SP or Sr. Consequently, the magnitude of either x or x would increase over the values obtained with the normal data.

The squared prediction error (SPE) indicates the difference between a sample and its projection into the k components retained in the model. Mathematically, the estimation error is given by

E = x - x = x = (I - PP T ) x (9)

And its energy as

SPE = ||x||2 = ||(I - PPT)x||2 . (10)

The process is considered normal if

SPE <S2 (11)

Where S2 is a confidence limit for SPE.

4.2 Wavelet analysis

The wavelet Transform is defined as the integral of the signal x(t) multiplied by a scaled and shifted version of a basic wavelet function \p(t), which satisfies the admissibility criteria [6],

c(a,b) = J s(t)^= —b)dt a e R+ - {o} b e R , (13)

JR Va a ,

where a is the so-called scaling parameter, b is the time localization parameter. Both a and b can be continuous or discrete variables. Multiplying each coefficient by an appropriately scaled and shifted wavelet yields the constituent wavelets of the original signal. For signals of finite energy, continuous wavelet synthesis provides the reconstruction formula,

s(t) =— f fc(a,b)-1 ^(t-b)db (14) Kw RR+ 4a b a

associated with the wavelet, which is used to define the details (high scale = low frequency content) in the decomposition, a scaling function \jdt), is used to define the approximation (low scale = high frequency content). Note that J^( x) = 1 while JV( x)dx = 0.To avoid intractable

computations when operating at every scale of the CWT, scales and positions can be chosen on a power of two, i.e. dyadic scales and positions. This defines the discrete wavelet transform (DWT). In this scheme, a and bare given by: (j,k) e Z2 : a = 2j,

b = k2J, Z = {o, ± 1, ±2, •••} (15)

Let us define:

(J, k) e Z2 : j = 2-j/2^(2-jt -k)

j (t) = 2- J /2ф(2- jt - k), (16)

A wavelet filter with impulse response g, which plays the role of the wavelet ^ , and a scaling filter with impulse response h, which plays the role of scaling function ф. Filters g and h are defined on a regular grid AZ where A is the sampling period (A=1). Then the discrete wavelet analysis can be described mathematically as

C (йу b) = c( j, k) = £ s(n) gj k(n)

neZ

a =2j, b=k2j, jeN,keN (17)

and the discrete wavelet synthesis as

s(t) = EEc(j,k)Wj,k(t). (18)

jeZ keZ

The detail at level j is defined as

R&I, 2012, №4

33

D(t) = X c( j,k )Wjk (t),

(19)

and the approximation at level j, Aj-1 = XDj . Obviously,

j>j

the following equations hold:

A j-1 = Aj + Dj

s = Aj +X Dj J

In this research work, Daubechies-6 wavelet is used for vibration signal processing and analysis.

(20)

4.3 MSPCA formulation

The recent success of wavelets and multi-scale methods in analysing and diagnosing multivariate processes suggest similar successes in the investigation of fault detection methods. In this work, multi-scale principal component analysis (MSPCA) is used for fault detection and diagnosis. MSPCA simultaneously extracts both the cross correlation across the sensors (PCA approach) and the auto-correlation within a sensor (wavelet approach).Using wavelets, the individual sensor signals are decomposed into approximations and details at different scales. Contributions from each scale are collected in separate matrices, and a PCA model is then constructed to extract correlation at each scale. The multi-scale nature of MSPCA formulation makes it suitable to work with process data that are typically non-

stationary and represent the cumulative effect of many underlying process phenomena, each operating at a different scale.

For the MSPCA formulation, consider a n x m data matrix X having m variables and n samples. Each of the m columns are first decomposed individually by applying a discrete wavelet transform (DWT). It may be noted that the same wavelet transform with the same level of decomposition L, is applied to each of the m variables. The wavelet approximation AL from each of the m decompositions is collected in one matrix of size m x n/ 2L (as shown by the AL solid line matrix in Fig.4).Similarly, the wavelet details (D1 to DL from each of the L levels) from each of the decompositions are collected in L corresponding matrices (as shown by the DL dashed and D1 dotted line matrices in Fig.4), with matrix size varying m x n/ 2i, i = 1,2, l , l . Thus a total of L+1matrices are formed, each being represented at a different scale, and captured trends within the sensors at the corresponding scale.

The PCA is then applied to each of the L+1matrices, the objective being to extract the correlation across the sensors.

Fig. 4. Multi-scale PCA,______wavelet approximations mode,.......wavelet details models

V. Fault detection with the Spectrum Energy (SE)

To compute the spectrum energy, a moving data window goes through the vibration wavelet detail coefficients shifting at a time according this expression:

nw

Ew =X[Sw(k)]2,

k=1

where Sw (k) is the k th spectrum coefficient within the

wth window and nw is the window length. To do this detection, the spectrum energy(SE) is applied on the detail

coefficients at different levels and comparing the level of this energy between the healthy and the defect bearing.

VI. Results and discussions

Experiments have been conducted on several signals of healthy and defective bearings with fault severity of different levels. However, due to limited space, results concerned with inner race faults only are considered. The acquired signals of the rotating machine in the normal operating mode and defective mode are shown in Fig.5 and 6. The simulation is done using two methods, the SPE based on the MSPCA and the spectrum Energy (SE).

34

R&I, 2012, №4

The detail coefficients (D1-D4) of the healthy bearing and the inner race defects are illustrated in Fig.7, 8. With regard to the MSPCA method, it is easily seen that the presence of failure provides a very important change on the detail coefficients (D1-D4) between the healthy bearing and the defect one. The presence of failure provides also a change in the correlation between the variables indicating an abnormal situation.

In this case, the projection of the measurement vector in the residual subspace will be increased. To detect this change, the squared Prediction Error (SPE) is used.

The SPE is compared with the limit corresponding to 95% of maximal amplitude of the healthy signal in Fig. 9.

According to figures (Figs. 9-16),it is seen that the magnitude of the eigenvalue is increased in the defective case relatively with the healthy operating mode, and the defective bearing showsan SPE 95% above the limit, which is strong indication of the fault presence.

Then to detect the failures with the SE, we have calculatedthe spectrum energy of the detail coefficients (D1-D4) and the results are shown in figures (Figs. 17-21). These results demonstrate that in the healthy mode operating (Fig. 17) the spectrum energy of the D3 is the most important magnitude relatively to the other coefficients. However, in the defective mode, the fault is localized in the high frequencies and the spectrum energy of the D1 is the most important magnitude in relation to the other coefficients (Figs. 18-21).

In conclusion, it can be said that the SPE can be used as an operational status indicator to discriminate between a safe operational mode and a defective one.

Vibration Signals of Healthy Bearing

500 1000 1500 2000 2500

0 Ь I

-0.5 C-------------1------------L-------------1-------------1- '

500 1000 1500 2000 2500

0.5

-0.5

0 c

_l_

_l_

_l_

500 1000 1500 2000 2500

3000

0 500 1000 1500 2000 2500 3000

Samples

Fig. 5. Output sensors: (a, b, c, d) of healthy Bearing

Fig. 6. Output sensors: (a, b, c, d) of defective Bearing

Fig. 7. Details coefficients (D1-D4)of healthy bearing

Fig. 8. Details coefficients (D1-D4) of defective bearing

R&I, 2012, №4

35

SPE(APd5) SPE(AP5) M » SPE(APd5) ш SPE(AP5)

Fig. 9 SPE Evolution of D1: of healthy, b - of defect bearings

Fig. 10 Eigenvalue evolution of D1: a - of healthy, b - of defect bearings

:PE evolution of Detail2 of healthy bearings Eigen value evolution of healthy bearing

100 200 300 400 500 1234

Time (ms) (a) (a)

5PE evolution of Detail2 of defect bearings Eigen value evolution of defect bearing

Fig.11. SPE Evolution of D2: a - of healthy, b - of defective bearings

Fig.12. Eigenvalue evolution of D2: a - of healthy, b - of defective bearings

36

R&I, 2012, №4

SPE evolution of Detail3 of healthy bearings

Time [ms) [a)

SPE evolution of Detail3 of defect bearings

Time [ms) [b)

Fig.13. SPE Evolution of D3: a - of healthy, b - of defective bearings

Eigen value evolution of healthy bearing

12 3 4

[a)

Eigen value evolution of defect bearing

12 3 4

(b>

Fig.14. Eigenvalue evolution of D3: a - of healthy, b - of defective bearings

Fig. 15. SPE Evolution of D4: a - of healthy, b - of defective bearings

Fig. 16. Eigenvalue evolution of D4: a - of healthy, b - of defective bearings

R&I, 2012, №4

37

Fig. 18. Energy of detail coefficient (D1-D4) of defective bearing (sensor n°1)

Fig. 19. Energy of detail coefficient (D1-D4) of defective bearing (sensor n°2)

Fig. 20. Energy of detail coefficient (D1-D4) of defective bearing (sensor n°3)

Fig. 21. Energy of detail coefficient (D1-D4) of defective bearing (sensor n° 4)

VII. Conclusion

In this work, the multi-scale principal component and spectral energy methods have been presented. To detect the presence of faults, the squared Prediction Error (SPE) and the spectrum energy of the detail coefficients are used. The SPE is compared with the limit corresponding to 95% of maximum amplitude of the healthy signal and the spectrum energy calculated from the detail coefficients. If the SPE is below the 95% limit, the sample is assumed to be in a normal operation, but if it is above, it indicates an abnormal situation. The results obtained from the SPE criteria based on the MSPCA method were confirmed by the Spectral Energy method. Then we can easily use the detail coefficient at level 1 as indicator of presence of the fault.

In fact, fault detection can be achieved by the multi-scale statistic SPE. From this technique, a considerable amount of information is acquired which can be used for the identification of different faults existing in rotating machines.

38

R&I, 2012, №4

In conclusion, it can be said that SPE can be used as an operational status indicator to discriminate between a safe operational mode and a defective one.

References

[1] Lachouri A.; Baiche K.; Djeghader R.; Doghmane N.; Ouchtati S. 2008.Analyze and fault diagnosis by multi-scale PCA. ICTTA, Damascus, page(s):1-6 ISBN: 978-1-4244-1751-3.

[2] Eren L.; DevaneyM. J.2004. Bearing damage detection via wavelet packet decomposition of the stator current, IEEE Trans on Instrumentation and Measurement, vol. 53, N. 2.

[3] Vasylius M.; Didziokas R.; Mazeika P.;Barzdaitis V.2008.The rotating system vibration and diagnostics, Mechanika-Kaunas: Technologija, 2008, Nr.4(72), p.54-58.

[4] Antoni J.; Randal R. B. 2006. The spectral kurtosis: application to the vibratory surveillance and diagnostics of rotating machines. Elsevier. Mechanical Systems and Signal Processing Vol.20,308-331.

[5] Castejon C.; Lara O.; Garcfa-Prada J. C. 2010. Automated diagnosis of rolling bearings using MRA and neural networks. Science direct, Elsevier. Mechanical Systems and Signal Processing, Vol.10.

[6] Baiche K.; Lachouri A.;Zelmat M. 2006. Vibration analysis of rotating machine by high resolution methods. Journal of Engineering and Applied Sciences, 1(4), pages 335-340.

[7] HarkatM. F. 2003.Faults detection and localisation by principal component analysis”. Doctorat thesis, Automatic Researsh Center of Nancy- UMR-CNRS 7039.

[8] Wei-Xin Ren.; Zeng-Shou Sun. 2008. Structural damage identification by using wavelet entropy, Engineering Structures. Elsevier, Engineering Structures, Vol. 30.

[9] Antoni J. 2006. The spectral kurtosis: a useful tool for characterising non-stationary signals. Elsevier. Mechanical Systems and Signal Processing 20, 282-307.

[10] Manish Misra; Henry YueH.; Joe Qin S.; Cheng Ling. 2002. Multivariate process monitoring and fault diagnosis by multi-scale PCA. Computers and Chemical Engineering 26, 1281-1293.

[11] Santoso S.; PowersE. J.; GradyW. V.; ParsonsA. C. 2000. Power quality disturbance waveform recognition using wavelet based neural classifier-Part 1: Theoretical foundation”. IEEE Trans Power Delivery, vol 15, pp. 222-228.

[12] Castro O. J. L.; SisamonC. C.;PradaJ.C.G. 2006. Bearing fault diagnosis based on neural network classification and wavelet transform.WAMUS'06 Proceedings of the 6th WSEAS international conference on Wavelet analysis & multirate systems.

R&I, 2012, №4

39