# A Novel Method for Mechanical Fault Diagnosis Based on Variational Mode Decomposition and Multikernel Support Vector Machine

Authors: Zhongliang Lv; Baoping Tang; Yi Zhou; Chuande Zhou
Journal: Shock and Vibration
 (2016)

---

## Abstract

A novel fault diagnosis method based on variational mode decomposition (VMD) and multikernel support vector machine (MKSVM) optimized by Immune Genetic Algorithm (IGA) is proposed to accurately and adaptively diagnose mechanical faults. First, mechanical fault vibration signals are decomposed into multiple Intrinsic Mode Functions (IMFs) by VMD. Then the features in time-frequency domain are extracted from IMFs to construct the feature sets of mixed domain. Next, Semisupervised Locally Linear Embedding (SS-LLE) is adopted for fusion and dimension reduction. The feature sets with reduced dimension are inputted to the IGA optimized MKSVM for failure mode identification. Theoretical analysis demonstrates that MKSVM can approximate any multivariable function. The global optimal parameter vector of MKSVM can be rapidly identified by IGA parameter optimization. The experiments of mechanical faults show that, compared to traditional fault diagnosis models, the proposed method significantly increases the diagnosis accuracy of mechanical faults and enhances the generalization of its application.


---

## Body

## 1. Introduction

To ensure the safe and reliable operation of mechanical equipments, vibrations are usually analyzed to diagnose mechanical faults [1–3]. Accurate diagnosis helps to make reasonable maintenance decision. However, the diagnosis accuracy is generally low and manual intervention of diagnosis is usually needed [4, 5], due to the features of large rotating machinery, including multicomponent coupling vibration, strong vibration noise interference, and the instability and nonlinearity of signals, as well as the low noticeability of early fault signals.As an adaptive method which processes signal in time-frequency domain, empirical mode decomposition (EMD) [6] can decompose the complicated and unstable signal into several nearly stable IMFs. Therefore, it has been widely applied to the diagnosis of mechanical fault [7]. For example, Loutridis [8] used EMD method for gear fault diagnosis. Cheng et al. used EMD method to diagnose bearing failures [9]. However, EMD method is essentially a binary filter bank. The frequency domain splitting feature of EMD makes it disadvantageous in dealing with fault signals [10]. Since the band center and bandwidth of fault signal are unknown, strong interference may be introduced if the fault signal falls in the broad bands of the first component. If the fault signal is in the bands of higher-order components, important information on the feature may be missed as the signal is filtered out by the narrow band of low-order components. Dragomiretskiy and Zosso, in 2014, proposed a new adaptive signal processing method called VMD [11]. This method is able to determine the frequency center and bandwidth of each component in the process of acquiring decomposed components by iteratively searching the optimal solution of variational models, thus adaptively realizing the frequency domain split and the effective separation of each component.In pattern recognition, support vector machine (SVM) [12] is based on Vapnik-Chervonenkis dimension theory and structural risk minimization principle. SVM finds the optimal compromise between model complexity and learning ability using limited sample information. It overcomes the drawback that traditional machine learning models are easy to get trapped in local minima. It has enormous potential to accurately classify the faults into multiple levels. Therefore, SVM has been widely applied to all kinds of nonlinear pattern recognition problems [13]. However, in complicated cases, especially when the data are heterogeneous [14] and samples are unevenly distributed [15], or samples are in large scales [16], SVM [17] begins to lose its advantages in accomplishing the tasks. MKSVM is a new machine learning model which combines all individual kernels by weights based on traditional single kernel SVM. MKSVM inherits the generalization and learning ability of single kernel SVM. Meanwhile, it reasonably adjusts the weight of each individual kernel and improves the adaptability and robustness of single kernels [18]. However, in the fault identification process using MKSVM, the identification performance of MKSVM is directly influenced by the choice of the function parameters of individual kernels and their weights. Traditional methods, including trial and error or traversing optimization, are not only complicated in computation, but also unable to acquire global optimal solution. Therefore, the adaptive diagnosis ability of MKSVM still needs further improvements. IGA [19] finds optimal solutions by synthetically considering the information interaction between antibodies of populations. Based on genetic algorithm [20], it integrates a series of mechanisms of biological immune system, such as antigen recognition, antibody diversity, density control, and elitist strategy. It greatly helps to avoid immature convergences and meanwhile preserves the global stochastic parallel searching character of genetic algorithm at the same time [19, 21]. Taking computational efficiency, stability, and global optimality factors into consideration, IGA is used to optimize the penalty parameter, weight factor, and kernel parameters of MKSVM, in order to improve the accuracy and stability of fault diagnosis, as well as to enhance the applicability of MKSVM.The rest of this paper is organized as follows. In Section2, the theories of EMD and VMD are introduced. In Section 3, the theories of SVM and MKSVM are firstly reviewed. Then, the IGA optimized MKSVM (IGA-MKSVM) method is presented. In Section 4, the fault diagnosis strategy based on VMD and IGA-MKSVM is discussed. Finally, the effectiveness of the proposed method is verified by experiments.


## 2. VMD Method

### 2.1. EMD

IMF components obtained by EMD should meet the following criteria: (1) In a data sequence, the number of extreme points and the number of zero crossing points are equal, or up to a difference of 1; (2) at any data point, the average of the local maximum envelope and the local minimum envelope is 0. The basic algorithm of EMD is as follows.Step 1.
First, determine the local extremes of signalx(t). Then use a cubic spline to connect all the local maximum points to form the upper envelope. After that, use another cubic spline to connect all the local minimum points to form the lower envelope. The upper and lower envelopes should enclose all data points. The average of the upper and lower envelopes is denoted as m(t). h(t) refers to the new signal, which is obtained by subtracting m(t) from x(t):(1)ht=xt-mt.
Repeat Step1 by K times until h1k(t) becomes the basic IMF component.Step 2.
Definec1(t), c1(t)=h1k(t), which is the first mode component obtained by processing the raw data. It should contain the shortest periodic component of raw signals. The residual component r1(t) can be calculated as follows:(2)r1t=xt-c1t.Step 3.
Since the residual componentr1(t) still contains information of long periodic components, r1(t) is still treated as new signal data. Repeat the above steps for the residual component rj(t), and the following results can be obtained:(3)r1t-c2t=r2t,r2t-c3t=r3t,⋮rN-1t-cNt=rNt.Step 4.
The original signalx(t) is finally decomposed into the sum of several IMF components ci(t) and a residual rN(t):(4)xt=∑i=1Ncit+rNt.


### 2.2. VMD

In VMD algorithm [11, 22], an intrinsic mode function is redefined as an AM-FM signal, which is expressed by(5)ukt=Aktcosφkt,where Ak(t) is the instantaneous amplitude of uk(t); wk(t) is the instantaneous frequency of uk(t) and wk(t)=φk′(t)=dφ(t)/dt. Ak(t) and ωk(t) are slowly varying with respect to the phase ϕk(t). That is, in the interval of [t-δ,t+δ] (where δ=2π/ϕk′(t)), uk(t) can be viewed as a harmonic signal with amplitude Ak(t) and frequency ωk(t).To obtain IMF components, VMD algorithm does not use the cycled screening stripping signal processing mode of EMD. Instead, VMD moves the signal decomposition process into the variational framework. It realizes adaptive signal decomposition by searching the optimal solution of the constrained variational model. The frequency center and bandwidth of each IMF component are updated in the iterative solving process of the variation model. The signal band is adaptively split according to the frequency domain features of the signal. Finally the narrow band IMF components are obtained.Assuming the original signalx(t) is decomposed into K IMF components, the corresponding constrained variational model is expressed as follows:(6)minuk,wk∑k∂tσt+jπt∗ukte-jwkt22s.t.∑kuk=xt,where uk=u1,…,uk represents the K IMF components decomposed by VMD; wk=w1,…,wk represents the frequency centers of all IMF components.To obtain the optimal solution of the above constrained variational problem, the following augmented Lagrange function is introduced:(7)Luk,wk,λ=α∑k∂tσt+jπt∗ukte-jwkt22+xt-∑kukt22+λt,xt-∑kukt,where α is a penalty factor, and λ is the Lagrange multiplier.The optimal solution of the constrained variational model is derived by using alternating direction multiplier algorithm, which solves the saddle point of the above augmented Lagrange function. The decomposition of the original signalx(t) is then obtained. The detailed implementation steps are as follows:(1)
Initializeuk1, wk1, λ1, and n as 0.(2)
n=n+1; perform the entire cycle.(3)
Perform the first inner cycle and updateuk according to:(8)ukn+1=argminukLui<kn+1,ui≥kn,win,λn.(4)
k=k+1; repeat step (3) until k=K and finish the first inner cycle.(5)
Perform the second inner cycle and updatewk according to:(9)wkn+1=argminLuin+1,wi<kn+1,wi<kn+1,wi>kn,λn.(6)
k=k+1; repeat step (5) until k=K. Finish the second inner cycle.(7)
Updateλ according to λn+1=λn+τ(f-∑ukn+1).(8)
Repeat steps (2) to (7) until the criterion∑ukn+1-ukn22/∑ukn22<ε is satisfied. Finish the entire cycle and output the results of K narrow band IMF components.


## 2.1. EMD

IMF components obtained by EMD should meet the following criteria: (1) In a data sequence, the number of extreme points and the number of zero crossing points are equal, or up to a difference of 1; (2) at any data point, the average of the local maximum envelope and the local minimum envelope is 0. The basic algorithm of EMD is as follows.Step 1.
First, determine the local extremes of signalx(t). Then use a cubic spline to connect all the local maximum points to form the upper envelope. After that, use another cubic spline to connect all the local minimum points to form the lower envelope. The upper and lower envelopes should enclose all data points. The average of the upper and lower envelopes is denoted as m(t). h(t) refers to the new signal, which is obtained by subtracting m(t) from x(t):(1)ht=xt-mt.
Repeat Step1 by K times until h1k(t) becomes the basic IMF component.Step 2.
Definec1(t), c1(t)=h1k(t), which is the first mode component obtained by processing the raw data. It should contain the shortest periodic component of raw signals. The residual component r1(t) can be calculated as follows:(2)r1t=xt-c1t.Step 3.
Since the residual componentr1(t) still contains information of long periodic components, r1(t) is still treated as new signal data. Repeat the above steps for the residual component rj(t), and the following results can be obtained:(3)r1t-c2t=r2t,r2t-c3t=r3t,⋮rN-1t-cNt=rNt.Step 4.
The original signalx(t) is finally decomposed into the sum of several IMF components ci(t) and a residual rN(t):(4)xt=∑i=1Ncit+rNt.


## 2.2. VMD

In VMD algorithm [11, 22], an intrinsic mode function is redefined as an AM-FM signal, which is expressed by(5)ukt=Aktcosφkt,where Ak(t) is the instantaneous amplitude of uk(t); wk(t) is the instantaneous frequency of uk(t) and wk(t)=φk′(t)=dφ(t)/dt. Ak(t) and ωk(t) are slowly varying with respect to the phase ϕk(t). That is, in the interval of [t-δ,t+δ] (where δ=2π/ϕk′(t)), uk(t) can be viewed as a harmonic signal with amplitude Ak(t) and frequency ωk(t).To obtain IMF components, VMD algorithm does not use the cycled screening stripping signal processing mode of EMD. Instead, VMD moves the signal decomposition process into the variational framework. It realizes adaptive signal decomposition by searching the optimal solution of the constrained variational model. The frequency center and bandwidth of each IMF component are updated in the iterative solving process of the variation model. The signal band is adaptively split according to the frequency domain features of the signal. Finally the narrow band IMF components are obtained.Assuming the original signalx(t) is decomposed into K IMF components, the corresponding constrained variational model is expressed as follows:(6)minuk,wk∑k∂tσt+jπt∗ukte-jwkt22s.t.∑kuk=xt,where uk=u1,…,uk represents the K IMF components decomposed by VMD; wk=w1,…,wk represents the frequency centers of all IMF components.To obtain the optimal solution of the above constrained variational problem, the following augmented Lagrange function is introduced:(7)Luk,wk,λ=α∑k∂tσt+jπt∗ukte-jwkt22+xt-∑kukt22+λt,xt-∑kukt,where α is a penalty factor, and λ is the Lagrange multiplier.The optimal solution of the constrained variational model is derived by using alternating direction multiplier algorithm, which solves the saddle point of the above augmented Lagrange function. The decomposition of the original signalx(t) is then obtained. The detailed implementation steps are as follows:(1)
Initializeuk1, wk1, λ1, and n as 0.(2)
n=n+1; perform the entire cycle.(3)
Perform the first inner cycle and updateuk according to:(8)ukn+1=argminukLui<kn+1,ui≥kn,win,λn.(4)
k=k+1; repeat step (3) until k=K and finish the first inner cycle.(5)
Perform the second inner cycle and updatewk according to:(9)wkn+1=argminLuin+1,wi<kn+1,wi<kn+1,wi>kn,λn.(6)
k=k+1; repeat step (5) until k=K. Finish the second inner cycle.(7)
Updateλ according to λn+1=λn+τ(f-∑ukn+1).(8)
Repeat steps (2) to (7) until the criterion∑ukn+1-ukn22/∑ukn22<ε is satisfied. Finish the entire cycle and output the results of K narrow band IMF components.


## 3. IGA-MKSVM

### 3.1. SVM

Assume two linearly separable sample sets,{(xi,yi),i=1,…,n},xi∈Rn,yi∈{+1,-1}. The general form of decision functions is fx=ω·x+b, and the decision surface equation is as follows: (10)ω·x+b=0.Normalizing the decision equation so that the samples of both classes satisfyf(x)≥1. The samples closest to the decision surface satisfy f(x)=1. The decision surface correctly classifies all samples; that is,(11)yiω·xi+b-1≥0,i=1,…,n.Support vectors are samples satisfying (11) and such that (1/2)ω2 is minimum. Those samples are on the lines of H1 and H2, as shown in Figure 1, where H is the optimal classification surface.Figure 1
Schematic of the optimal surface.Under the constraint of (11), the problem of optimal classification surface can be expressed as follows:(12)s.t.yiω·xi+b-1≥0,i=1,…,n,ϕω=12ω2.An optimal hyper surface is found such that the average classification error for the entire training sample set reaches minimum. Introducing a nonnegative relaxation factorξ and allowing the existence of misclassified samples, the decision surface ω·x+b=0 satisfies(13)yiω·x+b≥1-ξi,i=1,…,n.A penalty termc∑i=1nξi is then added to (1/2)ω2, and the following objective function is introduced:(14)ϕω,ξ=12ω2+c∑i=1nξi,where c is the penalty factor.Using Lagrange optimization method, the above optimal classification surface problem is converted into the following dual problem of convex quadratic programming optimization:(15)max∑i=1nai-12∑i=1n∑j=1naiajyiyjxi·xjs.t.ai≥0,i=1,…,n∑i=1naiyi=0.By solving the above problem, the optimal decision function can be described as follows:(16)fx=sgnω·x+b=sgn∑i=1nai∗yixi·x+b∗.


### 3.2. MKSVM

Since different kernel functions correspond to different decision functions, the selection of kernel functions is very important in fault identification using SVM, and it directly affects the identification accuracy of SVM.The kernel functions of SVM mainly include local kernel functions and global kernel functions. The Gaussian kernel function is a typical kernel function, which is described as follows:(17)KRBFxi,xj=exp-xi-xj2σ2.Polynomial function is a typical global kernel function, which is described as follows:(18)Kployxi,xj=xiT-xj+1d.Local kernel functions have strong learning ability but weak generalization ability; while global kernel functions have strong generalization ability but weak learning ability. In order to achieve better learning and generalization abilities of SVM, MKSVM is constructed based on local kernels and global kernels:(19)Kmixxi,xj=λKRBFxi,xj+1-λKployxi,xj,where xi and xj are feature vectors of input space; σ and d are kernel parameters of Gaussian kernel and polynomial kernel, respectively; λ(0<λ<1) is a weight factor. MKSVM combines the advantages of each single kernel and has better classification and identification performance.


### 3.3. IGA-MKSVM

IGA treats the object problem to be solved as biological invasion antigen and the feasible solution of the problem as antibody. The searching process of the optimal solution can be viewed as the process of seeking maximum antigen affinity antibodies by biological systems. The inhibition and promotion of antibodies can ensure the diversity of antibodies in the population and improve the local searching ability of GA. Crossover and mutation of antibodies can ensure that the antibody population evolves towards the direction of high fitness and maintain the diversity of the population. The memory unit accelerates searching by constantly updating with better solutions, which improves the global searching capability of the algorithm. The flowchart of IGA algorithm is shown in Figure2. Refer to [18, 19] for more details of IGA.Figure 2
Flowchart of IGA.In IGA-MKSVM, IGA algorithm is used to optimize the weight factor, penalty parameter, and kernel parameter. First, an antibody gene vectorX is defined, which is composed of the weight factor of MKSVM λ, penalty parameter c, and parameters σ and d of kernel function:(20)X=λ,c,σ,d.To minimize the square error between the actual output and the expected output of MKSVM, the fitness functionf(xi) of MKSVM is defined based on the classification accuracy of training sample E(xi):(21)fxi=Exi.The flowchart of IGA-MKSVM algorithm is shown in Figure2, and the detailed steps are as follows:(1)
Initialize the population and determine the population size, fitness threshold, and maximum iteration number. Determine the initial vector of each antibody within the ranges of penalty parameter and kernel parameters.(2)
Compute the fitness value of each antibody according to (21) using the current antibody and training samples.(3)
For the current population, choose the antibody with the highest fitness value as the elite antibody, and save this antibody in a special variable.(4)
If it is the first generation of antibody population, go to step (7); otherwise, go to the next step.(5)
Determine the fitness value of each antibody. If no antibody in the current antibody generation has the same fitness value as the elite antibody, then replace the antibody having the smallest fitness value in the current antibody population by the elite antibody saved in a special variable; otherwise, go to the next step.(6)
If the maximum fitness value in the current antibody population is larger than that of the elite antibody, then copy the antibody with the maximum fitness value into the special variable to replace the current elite antibody; otherwise, go to the next step.(7)
According to the similarity definition, compute the density and selection probability of each antibody; perform selection and copy operations for the antibody population according to the selection probability.(8)
Perform crossover and mutation operation for the antibody population.(9)
Judge whether the criterion of termination is satisfied. If yes, output the results and the algorithm ends; otherwise, return to step (2) and continue the cycle.


## 3.1. SVM

Assume two linearly separable sample sets,{(xi,yi),i=1,…,n},xi∈Rn,yi∈{+1,-1}. The general form of decision functions is fx=ω·x+b, and the decision surface equation is as follows: (10)ω·x+b=0.Normalizing the decision equation so that the samples of both classes satisfyf(x)≥1. The samples closest to the decision surface satisfy f(x)=1. The decision surface correctly classifies all samples; that is,(11)yiω·xi+b-1≥0,i=1,…,n.Support vectors are samples satisfying (11) and such that (1/2)ω2 is minimum. Those samples are on the lines of H1 and H2, as shown in Figure 1, where H is the optimal classification surface.Figure 1
Schematic of the optimal surface.Under the constraint of (11), the problem of optimal classification surface can be expressed as follows:(12)s.t.yiω·xi+b-1≥0,i=1,…,n,ϕω=12ω2.An optimal hyper surface is found such that the average classification error for the entire training sample set reaches minimum. Introducing a nonnegative relaxation factorξ and allowing the existence of misclassified samples, the decision surface ω·x+b=0 satisfies(13)yiω·x+b≥1-ξi,i=1,…,n.A penalty termc∑i=1nξi is then added to (1/2)ω2, and the following objective function is introduced:(14)ϕω,ξ=12ω2+c∑i=1nξi,where c is the penalty factor.Using Lagrange optimization method, the above optimal classification surface problem is converted into the following dual problem of convex quadratic programming optimization:(15)max∑i=1nai-12∑i=1n∑j=1naiajyiyjxi·xjs.t.ai≥0,i=1,…,n∑i=1naiyi=0.By solving the above problem, the optimal decision function can be described as follows:(16)fx=sgnω·x+b=sgn∑i=1nai∗yixi·x+b∗.


## 3.2. MKSVM

Since different kernel functions correspond to different decision functions, the selection of kernel functions is very important in fault identification using SVM, and it directly affects the identification accuracy of SVM.The kernel functions of SVM mainly include local kernel functions and global kernel functions. The Gaussian kernel function is a typical kernel function, which is described as follows:(17)KRBFxi,xj=exp-xi-xj2σ2.Polynomial function is a typical global kernel function, which is described as follows:(18)Kployxi,xj=xiT-xj+1d.Local kernel functions have strong learning ability but weak generalization ability; while global kernel functions have strong generalization ability but weak learning ability. In order to achieve better learning and generalization abilities of SVM, MKSVM is constructed based on local kernels and global kernels:(19)Kmixxi,xj=λKRBFxi,xj+1-λKployxi,xj,where xi and xj are feature vectors of input space; σ and d are kernel parameters of Gaussian kernel and polynomial kernel, respectively; λ(0<λ<1) is a weight factor. MKSVM combines the advantages of each single kernel and has better classification and identification performance.


## 3.3. IGA-MKSVM

IGA treats the object problem to be solved as biological invasion antigen and the feasible solution of the problem as antibody. The searching process of the optimal solution can be viewed as the process of seeking maximum antigen affinity antibodies by biological systems. The inhibition and promotion of antibodies can ensure the diversity of antibodies in the population and improve the local searching ability of GA. Crossover and mutation of antibodies can ensure that the antibody population evolves towards the direction of high fitness and maintain the diversity of the population. The memory unit accelerates searching by constantly updating with better solutions, which improves the global searching capability of the algorithm. The flowchart of IGA algorithm is shown in Figure2. Refer to [18, 19] for more details of IGA.Figure 2
Flowchart of IGA.In IGA-MKSVM, IGA algorithm is used to optimize the weight factor, penalty parameter, and kernel parameter. First, an antibody gene vectorX is defined, which is composed of the weight factor of MKSVM λ, penalty parameter c, and parameters σ and d of kernel function:(20)X=λ,c,σ,d.To minimize the square error between the actual output and the expected output of MKSVM, the fitness functionf(xi) of MKSVM is defined based on the classification accuracy of training sample E(xi):(21)fxi=Exi.The flowchart of IGA-MKSVM algorithm is shown in Figure2, and the detailed steps are as follows:(1)
Initialize the population and determine the population size, fitness threshold, and maximum iteration number. Determine the initial vector of each antibody within the ranges of penalty parameter and kernel parameters.(2)
Compute the fitness value of each antibody according to (21) using the current antibody and training samples.(3)
For the current population, choose the antibody with the highest fitness value as the elite antibody, and save this antibody in a special variable.(4)
If it is the first generation of antibody population, go to step (7); otherwise, go to the next step.(5)
Determine the fitness value of each antibody. If no antibody in the current antibody generation has the same fitness value as the elite antibody, then replace the antibody having the smallest fitness value in the current antibody population by the elite antibody saved in a special variable; otherwise, go to the next step.(6)
If the maximum fitness value in the current antibody population is larger than that of the elite antibody, then copy the antibody with the maximum fitness value into the special variable to replace the current elite antibody; otherwise, go to the next step.(7)
According to the similarity definition, compute the density and selection probability of each antibody; perform selection and copy operations for the antibody population according to the selection probability.(8)
Perform crossover and mutation operation for the antibody population.(9)
Judge whether the criterion of termination is satisfied. If yes, output the results and the algorithm ends; otherwise, return to step (2) and continue the cycle.


## 4. Fault Diagnostic Strategy

The flowchart of fault diagnosis based on VMD and IGA-MKSVM is shown in Figure3.(1)
First, decompose the obtained fault vibration signals of rotating machinery by VMD, and generateL IMF components that contain most fault information.(2)
Extract 5 statistical features in time domain and 7 statistical features in frequency domain. These features constitute aL×12 combined domain feature set with high dimension. The indices of the 5 time domain and 7 frequency domain features are presented in Tables 1 and 2.(3)
Use the SS-LLE algorithm to merge the combined high-dimensional feature set and reduce its dimension and then input the feature set to IGA-MKSVM as feature vectors.(4)
Use training samples to train the IGA-MKSVM and obtain the optimal weight factorλ, penalty parameter c, and kernel function parameters σ and d. Then, use the optimal parameter vector to reconstruct the fault diagnosis model and identify failure modes of test samples. Finally, output the diagnostic results.Table 1
The dimensionless time domain characteristic index.
Number
Feature expression
1
X
-
=
1
N
∑
i
=
1
N
x
i
2
X
rms
=
1
N
∑
i
=
1
N
x
i
2
3
β
=
1
N
∑
i
=
1
N
x
i
4
4
σ
x
2
=
1
N
-
1
∑
i
=
1
N
x
i
-
X
-
2
5
X
r
=
1
N
∑
i
=
1
N
x
i
2Table 2
The characteristic parameters in frequency domain.
Number
Feature expression
6
p
1
=
∑
k
=
1
K
s
k
K
7
p
2
=
∑
k
=
1
K
s
k
-
p
1
2
K
-
1
8
p
3
=
∑
k
=
1
K
s
k
-
p
1
3
K
p
2
3
9
p
4
=
∑
k
=
1
K
f
k
4
s
k
∑
k
=
1
K
f
k
2
s
k
10
p
5
=
∑
k
=
1
K
f
k
s
k
∑
k
=
1
K
s
k
11
p
6
=
∑
k
=
1
K
f
k
-
p
5
2
s
k
K
12
p
7
=
∑
k
=
1
K
f
k
2
s
k
∑
k
=
1
K
s
k
∑
k
=
1
K
f
k
4
s
k
s
(
k
) is the frequency spectrum of signal x(n), k=1,2,3,…,K,  K is the number of spectrum lines, and fk is the frequency of the kth spectrum.Figure 3
The process of fault diagnosis based on VMD and IGA-MSVM.


## 5. Examples of Application

### 5.1. Experiment Set-Up and Signal Acquisition

The experimental test system is composed of a speed motor, a driving belt, a coupling, a test bearing, a magnetic brake, an acceleration sensor, and a signal record analyzer. The bearing model is N205. The experimental devices are shown in Figure4.Figure 4
Experimental devices.Three N205 bearings were used to simulate different damages, including the outer ring damage, rolling damage, and inner ring damage. Bearing damages were implemented by processing a slot in different parts of bearings with a laser cutting machine. The width of the slot is 0.3 mm, and the depth is 0.1 mm. Different bearing faults are shown in Figure5. Four modes of bearing vibration signals were collected. They are (I) normal, (II) inner ring damage, (III) outer ring damage, and (IV) rolling body damage. The main experimental parameters are tabulated in Table 3.Table 3
Experimental parameters of rolling bearing.
Experimental parameter
Value
Outer ring diameter
52 mm
Inner ring diameter
25 mm
Rolling ring diameter
7.5 mm
The number of rolling elements
12
Bearing contact angle
0
°
Driving speed
1000 r/min
Bearing type
N205
Sampling frequency
10 KHzFigure 5
Rolling bearing fault.ICP accelerometers and data acquisition equipment DP/INV306U were used to collect the vibration signals of bearings. The sampling frequency is 10 kHz. In the input terminal, a low pass filter was used for antialiasing. For each operating status, 120 sample groups were collected, with a total of 10000 points for each sample group. Thirty sample groups were randomly selected as training samples, and the other 90 sample groups were test samples. Figure6 shows the time domain waveforms of the four modes.Figure 6
The original signal waveforms in different status of rolling bearing.
(a)
Normal operating status
(b)
Bearing outer race fault
(c)
Bearing inner race fault
(d)
Bearing rolling element faultVMD decomposition was performed for the four modes of vibration signals. The decomposition results of VMD and EMD for normal vibration signals of bearings are shown in Figure7. The decomposition process of vibration signals in other operating statuses is similar. The values of the statistical features in 5 time domains and 7 frequency domains were extracted for each operating status. The first 8 IMF components containing the most fault information were used as the feature vector of the diagnostic model. SS-LLE was used for dimension reduction. A three-dimensional sensitive feature index chart after dimension reduction is shown in Figure 8.Figure 7
The decomposition results of VMD and EMD for a normal vibration signal.
(a)
VMD decomposition
(b)
EMD decompositionFigure 8
Comparison of dimension reduction after the decomposition of VMD and EMD.
(a)
(b)


### 5.2. Experimental Results and Analysis

The fault samples after dimension reduction were inputted into the IGA-MKSVM diagnosis model for classification recognition. The recognition results are shown in Figure9. It is clear in Figure 9 that the average classification accuracy by using VMD decomposition is 94.99%. It is higher than that of EMD decomposition which is 91.38%. Therefore, in the following experiments, VMD decomposition was used to construct the combined domain feature sets.Figure 9
Average accuracy of classification by two different decomposition methods.To further verify the advantage and stability of the IGA-MKSVM fault diagnosis model, the standard MKSVM, SVM with Gaussian kernel (GSVM), and SVM with polynomial kernel (PSVM) were choose for comparison. We used fold cross validation [23] with a total of 5 times of validation. The synthetic diagnostic results are shown in Table 4. It is clear in Table 4 that IGA-MKSVM diagnosis model achieves the highest fault diagnosis accuracy. The recognition accuracy of the standard MKSVM is higher than that of SVM with a single kernel. Since the MKSVM combines the advantages of Gaussian kernel and polynomial kernel, it has higher recognition accuracy than SVM with a single kernel. Meanwhile, by obtaining the optimal control parameters through IGA, MKSVM model integrates the prior knowledge of fault features and improves the recognition performance.Table 4
Comparison of classification accuracy with four different parameters.
Diagnostic model
Average percentage of correct recognition (%)
Normal
Outer ring fault
Inner ring fault
Rolling element fault
IGA-MSVM
100
96.67
93.33
90.00
MSVM
100
86.67
91.11
85.56
GSVM
100
81.11
84.44
74.44
PSVM
100
78.89
82.22
68.89To verify the robustness and generalization ability of IGA-MKSVM diagnostic model, stochastic noise signals withα=0.1,0.2,0.3 were added to the original signal. The original feature vector then becomes X′=[x1′,x2′,…,x8′], where xk′=[1+α·rand(1)]xk, and xk is the original time domain feature of IMF. The noisy feature vector was inserted into the IGA-MKSVM diagnostic model, standard MKSVM, GSVM, and PSVM. The results of them were compared by using 5-fold cross validation. The diagnostic results are shown in Table 5, which indicate that the IGA-MKSVM diagnostic model not only has relatively high diagnosis accuracy, but also has excellent noise immunity. When the extracted fault features are interfered at certain extent by noise, the diagnosis accuracy is still relatively high. In practical engineering applications, due to the interference of complicated environmental factors, it is hard to acquire sufficient fault information. Standard MKSVM and SVM with a single kernel are inadequate to engineering applications in such cases. However, IGA-MKSVM diagnostic model has noise immunity and good generalization ability, which is applicable to real fault diagnosis in engineering applications.Table 5
Comparison of the noise immunity by using four different models for fault diagnosis.
Diagnostic model
Average correct ratio of recognition (%)
α
=
0.1
α
=
0.2
α
=
0.3
IGA-MSVM
95.56
91.11
80.00
MSVM
88.89
83.33
71.11
GSVM
81.11
78.89
63.33
PSVM
77.78
68.89
48.89


## 5.1. Experiment Set-Up and Signal Acquisition

The experimental test system is composed of a speed motor, a driving belt, a coupling, a test bearing, a magnetic brake, an acceleration sensor, and a signal record analyzer. The bearing model is N205. The experimental devices are shown in Figure4.Figure 4
Experimental devices.Three N205 bearings were used to simulate different damages, including the outer ring damage, rolling damage, and inner ring damage. Bearing damages were implemented by processing a slot in different parts of bearings with a laser cutting machine. The width of the slot is 0.3 mm, and the depth is 0.1 mm. Different bearing faults are shown in Figure5. Four modes of bearing vibration signals were collected. They are (I) normal, (II) inner ring damage, (III) outer ring damage, and (IV) rolling body damage. The main experimental parameters are tabulated in Table 3.Table 3
Experimental parameters of rolling bearing.
Experimental parameter
Value
Outer ring diameter
52 mm
Inner ring diameter
25 mm
Rolling ring diameter
7.5 mm
The number of rolling elements
12
Bearing contact angle
0
°
Driving speed
1000 r/min
Bearing type
N205
Sampling frequency
10 KHzFigure 5
Rolling bearing fault.ICP accelerometers and data acquisition equipment DP/INV306U were used to collect the vibration signals of bearings. The sampling frequency is 10 kHz. In the input terminal, a low pass filter was used for antialiasing. For each operating status, 120 sample groups were collected, with a total of 10000 points for each sample group. Thirty sample groups were randomly selected as training samples, and the other 90 sample groups were test samples. Figure6 shows the time domain waveforms of the four modes.Figure 6
The original signal waveforms in different status of rolling bearing.
(a)
Normal operating status
(b)
Bearing outer race fault
(c)
Bearing inner race fault
(d)
Bearing rolling element faultVMD decomposition was performed for the four modes of vibration signals. The decomposition results of VMD and EMD for normal vibration signals of bearings are shown in Figure7. The decomposition process of vibration signals in other operating statuses is similar. The values of the statistical features in 5 time domains and 7 frequency domains were extracted for each operating status. The first 8 IMF components containing the most fault information were used as the feature vector of the diagnostic model. SS-LLE was used for dimension reduction. A three-dimensional sensitive feature index chart after dimension reduction is shown in Figure 8.Figure 7
The decomposition results of VMD and EMD for a normal vibration signal.
(a)
VMD decomposition
(b)
EMD decompositionFigure 8
Comparison of dimension reduction after the decomposition of VMD and EMD.
(a)
(b)


## 5.2. Experimental Results and Analysis

The fault samples after dimension reduction were inputted into the IGA-MKSVM diagnosis model for classification recognition. The recognition results are shown in Figure9. It is clear in Figure 9 that the average classification accuracy by using VMD decomposition is 94.99%. It is higher than that of EMD decomposition which is 91.38%. Therefore, in the following experiments, VMD decomposition was used to construct the combined domain feature sets.Figure 9
Average accuracy of classification by two different decomposition methods.To further verify the advantage and stability of the IGA-MKSVM fault diagnosis model, the standard MKSVM, SVM with Gaussian kernel (GSVM), and SVM with polynomial kernel (PSVM) were choose for comparison. We used fold cross validation [23] with a total of 5 times of validation. The synthetic diagnostic results are shown in Table 4. It is clear in Table 4 that IGA-MKSVM diagnosis model achieves the highest fault diagnosis accuracy. The recognition accuracy of the standard MKSVM is higher than that of SVM with a single kernel. Since the MKSVM combines the advantages of Gaussian kernel and polynomial kernel, it has higher recognition accuracy than SVM with a single kernel. Meanwhile, by obtaining the optimal control parameters through IGA, MKSVM model integrates the prior knowledge of fault features and improves the recognition performance.Table 4
Comparison of classification accuracy with four different parameters.
Diagnostic model
Average percentage of correct recognition (%)
Normal
Outer ring fault
Inner ring fault
Rolling element fault
IGA-MSVM
100
96.67
93.33
90.00
MSVM
100
86.67
91.11
85.56
GSVM
100
81.11
84.44
74.44
PSVM
100
78.89
82.22
68.89To verify the robustness and generalization ability of IGA-MKSVM diagnostic model, stochastic noise signals withα=0.1,0.2,0.3 were added to the original signal. The original feature vector then becomes X′=[x1′,x2′,…,x8′], where xk′=[1+α·rand(1)]xk, and xk is the original time domain feature of IMF. The noisy feature vector was inserted into the IGA-MKSVM diagnostic model, standard MKSVM, GSVM, and PSVM. The results of them were compared by using 5-fold cross validation. The diagnostic results are shown in Table 5, which indicate that the IGA-MKSVM diagnostic model not only has relatively high diagnosis accuracy, but also has excellent noise immunity. When the extracted fault features are interfered at certain extent by noise, the diagnosis accuracy is still relatively high. In practical engineering applications, due to the interference of complicated environmental factors, it is hard to acquire sufficient fault information. Standard MKSVM and SVM with a single kernel are inadequate to engineering applications in such cases. However, IGA-MKSVM diagnostic model has noise immunity and good generalization ability, which is applicable to real fault diagnosis in engineering applications.Table 5
Comparison of the noise immunity by using four different models for fault diagnosis.
Diagnostic model
Average correct ratio of recognition (%)
α
=
0.1
α
=
0.2
α
=
0.3
IGA-MSVM
95.56
91.11
80.00
MSVM
88.89
83.33
71.11
GSVM
81.11
78.89
63.33
PSVM
77.78
68.89
48.89


## 6. Conclusions

To achieve accurate and adaptive identification of mechanical faults, a new fault diagnosis method based on VMD and IGA optimized MSVM is proposed in this paper. IGA overcomes the immature convergence problem of traditional algorithms and is able to solve the initial parameter selection problem of MSVM. From this aspect, it makes MSVM more applicable, robust, and accurate. Experimental results demonstrate that, compared to traditional methods, the combination of VMD method and IGA-MSVM for fault diagnosis produces more repeatable results, with stronger generalization ability as well as anti-interference ability, providing a new effective method for mechanical fault diagnosis.



---

*Source: 3196465-2016-01-06.xml*