# Intelligent Roller Bearing Fault Diagnosis in Industrial Internet of Things

Authors: Ji Xu; Hong Zhou; Yanjun Fang
Journal: Wireless Communications and Mobile Computing
 (2022)

---

## Abstract

Advanced research studies on industrial Internet of things require effective feature extraction and accurate machinery health state evaluation. For roller bearing, a well-known mechanical component most extensively used in the industry, its running status directly affects the operation of the entire machinery and equipment. For intelligent fault diagnosis of roller bearing, the selection of the intrinsic mode function (IMF) modes in approaches of ensemble empirical mode decomposition (EEMD)/variational mode decomposition (VMD) becomes a tricky problem. To solve this problem, this study proposed an efficient scheme on roller bearing fault diagnosis that combines the refined composite multivariate multiscale sample entropy (RCMMSE) with different classifiers. Firstly, the synthetic noise signals are introduced to compare the effectiveness of the multiscale sample entropy (MSE) and the RCMMSE models. Secondly, the random noise signals are used to compare the performance of EEMD and VMD methods, where the envelope spectrum characteristics of fault signals are well described. Moreover, EEMD/VMD methods are utilized to decompose the roller bearing vibration signals into various modes to get the entropy values. Finally, the obtained RCMMSE is adopted as a feature vector and subsequently employed as an input of support vector machine, random forest, and probabilistic neural network models to conduct roller bearing fault identification. The extensive experimental results prove that this proposed scheme performs well and the classification accuracy of VMD-RCMMSE is higher than EEMD-RCMMSE.


---

## Body

## 1. Introduction

With the rapid development of Internet of things (IoT) [1, 2] and Industry 4.0 [3], there are increasingly massive real-time data from various types of mechanical equipment [4]. The availability of these data that contain abundant information about machine health has attracted more and more enterprises’ attention. It has been proved that large volume, high velocity, and diversity mechanical big data are the major properties of mechanical big data [5, 6]. Effective feature extraction from these data and accurate machinery health state evaluation with ever-accelerated updating of schemes have become hot research issues in the prognostic and health management systems in the era of industrial IoT [7].For roller bearing, a well-known mechanical component most extensively used in the industry, its operating status directly affects the operation of the entire machinery and equipment. Roller bearing failure is an important factor leading to mechanical equipment failure, so timely detection and fault diagnosis are of great significance, and analysing the vibration signal collected by the sensors to determine its failure is a commonly adopted scheme [8]. As the roller bearings are usually working in the vibration source environment, with the destruction of complex forms, the vibration signal is typically nonlinear and nonstationary. How to excavate the fault feature based on the vibration signals has been a research hotspot [9, 10].The widely used signal decomposing schemes include wavelet analysis [11, 12] and EMD [13]. However, the wavelet transformation is essentially nonadaptive due to the configuration of wavelet basis and decomposing layers. In EMD, the complicated signal can be self-adaptively decomposed into some IMFs with a residual component. Nevertheless, mode mixing is a stumbling block in EMD. To deal with this challenge, EEMD is proposed [14], in which a complex signal can be disintegrated into IMFs in terms of the local time-scale characteristic of the signal. In recent, EEMD has been extensively employed in fault detection [15].Using the EMD or EEMD schemes for fault feature extraction has been widely concerned. However, these two schemes have some disadvantages. The recursive mode decomposition in EMD/EEMD will propagate the envelope estimation error continuously. The signal contains no noise or intermittent signal, which leads to the decomposition of the mode mixing. Although the white noise scheme was used to suppress the mode mixing, the scheme needs to be carried out several hundreds of times of EMD/EEMD operation and will break out the signal composed of more than a real component. Moreover, EMD and EEMD cannot be separated correctly from the close frequencies, and it is quite challenging to select the suitable number of IMFs when applying EMD and EEMD.To address the abovementioned problems, a neoteric scheme of signal decomposition estimation called VMD was proposed [16–18]. The whole framework of the alternating direction multiplier scheme was resorted to consecutively refresh the modes and their centre frequencies and gently demodulate the modes into the corresponding base frequency bands. Ultimately, each mode and corresponding centre frequencies were abstracted together. In comparison with the recursive filter pattern of EMD and EEMD, VMD converts the signal into variational and non-recursive decomposition patterns, and it consists of adaptive Wiener filter groups in nature. VMD can isolate two pure harmonic signals with similar frequencies [16–18]. The IMFs obtained by EEMD denote the natural oscillatory pattern inserted into the signals, where the entropy values of every IMF were often abstracted as a feature to discover the properties of the vibration signals [19, 20].Entropy was one of the schemes to evaluate the time-series complexity. S. M. Pincus proposed the approximate entropy (AE) [21]. However, the length of the time series makes a critical influence on the performance of AE scheme. Hence, the value of AE is conformably lower than the intended one and fails relative coherence particularly when the data length is short. Based on AE, SE was proposed to deal with this challenge [22]. It has the advantages of short data, stable and low noise and interference capacity, and good consistency in the region of large parameter range. It can be noticed that the irregularity of time series can be only reflected by AE and SE individually. When the roller bearings fail, not only the frequencies but also the corresponding complexity of vibration signals has tremendous deviations. Accordingly, ME can be viewed as a property index for fault diagnosis [23]. Considering the MSE, the feature of the vibration signals can be abstracted under all kinds of conditions, and the eigenvectors were counted as the input of adaptive neuro-fuzzy inference system (ANFIS) for roller bearing fault recognition [22]. Moreover, the SE, to some extent, was undefined as no standard vectors, which were in accordance with one another. Undefined or imprecise SE leads to the degradation of authenticity of MSE algorithm. To overcome this drawback, the RCMMSE was revolved by Zhang et al. [24] to overcome these challenges. It is demonstrated that RCMMSE can not only ascend the precision of entropy assessment but also descend the probability of inducing undefined entropy.However, after EEMD and VMD, a necessary step is how to select the number of modes. For example, the number of patterns was usually determined by the correlation coefficient and mutual information scheme between each component and the original signals [18, 25]. In addition, MSE and RCMMSE can only get the single-channel signal. Therefore, it cannot precisely display the overall signal information. The multivariate multiscale sample entropy (MMSE) was a scheme that gets different time series and conducts various embedding aspects, delay time, and amplitude ranges of data channels in a strict way. Hence, it can directly analyse multichannel data [26]. Therefore, a scheme that relied on VMD and RCMMSE was proposed in this study to get the vibration feature of roller bearings.With the advent of computer techniques, a lot of fault recognition algorithms, such as SVM, RF, and PNN [27, 28], have been broadly applied to fault diagnosis. Therefore, the SVM, RF, and PNN were resorted to fulfill the fault recognition.As mentioned above, a scheme based on RCMMSEVMD and VMD is designed. At first, the vibration signals were disintegrated into series IMF/BL-IMF modes by the EEMD/VMD. Secondly, the RCMMSE was utilized to figure out the entropy values. Finally, the SVM/RF/PNN models were utilized to achieve the fault recognition.The rest of the study is organized as follows. It is shown that the review of VMD and RCMMSE schemes, respectively, is presented in Section2. The comparison of MSE/RCMMSE and EEMD/VMD is presented in Section 3. Section 4 gives the procedures of the proposed scheme, experimental data sources, and parameter selection. Section 5 provides the experimental validation. Conclusions are given in Section 6.


## 2. Theoretical Framework of VMD and RCMMSE Models

### 2.1. Basic Principle of Variational Mode Decomposition

The VMD process is divided into establishment and solution of variational constraint problems. Assuming there is a limited bandwidth in each mode, the variational problem is formulated ask mode functions ukt, which could minimize the estimated bandwidth. The summation of each mode is constrained to be corresponded to the input signal f. In particular, the structure of this problem could be divided into the following steps.


#### 2.1.1. Variational Problem Formulation

①
Hilbert transform: the analytic signal of every mode function is obtained, whose purpose is to get its spectrum:(1)δt+jπt∗ukt.②
The centre frequencye−jωkt of each modal signal is estimated. The spectrum of each mode is modified to be adaptive to its own baseband.(2)δt+jπt∗ukte−jωkt.③
TheL2 is derived, and it denotes the signal gradient, and the signal’s bandwidth is obtained. Therefore, constrained variational questions can be obtained as follows:(3)minuk,ωk∑k∂tδt+jπtμkt−jωkt2s.t.∑kμk=f,where μ is the decomposition of the k mode components, μ=μ1,μ2,…,μk, ω is the centre frequencies of each mode, and ω=ω1,ω2,…,ωk.


#### 2.1.2. The Solution of the Problem

①
The above constraint problem could be modified as a nonbinding problem by adding it to the penalty factor∂ and the Lagrangian multiplication operator λt, and here, the second penalty factor is used to guarantee the accuracy of the signal with the Gaussian noise, and the Lagrangian operator is adopted to ensure the constraint condition satisfied strictly, and the extended Lagrangian description is as follows:(4)Luk,ωk,λ≔∂∑k∂tδt+jπt∗utke−jwkt2+ft−∑kukt22+λt,ft−∑kukt.②
VMD used multiplication operator alternating direction method of multipliers (ADMMs) to solve the above variable problem. By alternately updating A, B, and C to seek the “saddle” of the Lagrangian expression:(5)ukn+1=argminuk∈ka∂tδt+jπt∗ukte−jwkt22+ft−∑iuit+λt222,where ωk is equal to ωkn+1, and ∑tuit is equal to ∑i≠kuitn+1.Equation (5) is transformed to another domain by frequency by the Parseval/Plancherel Fourier equidistant transformation:(6)u^kn+1=argminajω1+sgnω+ωk⋅u^kω+ωk22+f^ω−∑iu^iω+λ^ω222.The first item ofω with ω−ωk instead is as follows:(7)u^kn+1=argminajω−ωk1+sgnωu^kω22+f^ω−∑iu^iω+λ^ω222.Equation (7) is converted into a nonnegative integral form with frequency interval:(8)u^kn+1=argminu^k,uk∈X∫0∞4aω−ωk2u^kw2+2f^ω−∑iu^iω+λ^ω22dω.In this case, the solution of the quadratic optimization problem is as follows:(9)u^kn+1ω=f^ω−∑i≠ku^iω+λ^ω/21+2aω−ωk2.Based on the above process, the centre frequency value is converted to the frequency domain:(10)ωkn+1=argmin∫0∞ω−ωk2u^kω2dω.The centre frequency updated scheme is as follows:(11)ωkn+1=∫0∞ωu^kω2dω∫0∞u^kω2dω,where u^kn+1ω is equivalent to the current remaining amount of f^ω−∑i≠ku^iω Wiener filtering. A is the central gravity of the spectrum of the current mode function. For inverse Fourier transmission of u^kω, then the real part is ukt.The procedure of VMD algorithm is as follows:(1)
Initialize theu^k1, ωk1, λ^1, and n.(2)
Updateuk and ωk followed by (9) and (11).(3)
Updateλ, and(12)λ^n+1ω←λ^nω+τf^ω−∑ku^kn+1ω.(4)
For a given accuracye>0, if ∑ku^kn+1−u^kn22/u^kn22<e, stop iteration; otherwise, return to Step 2.


### 2.2. Basic Principle of Refined Composite Multivariate Multiscale Sample Entropy

#### 2.2.1. Sample Entropy

The SE uses the measurement of the exponential function for two sequences with a tolerancer from m points to remain r of each other at the next point. For a time series with N sample points, the procedure of SE calculation is shown as follows:(1)
Them length vectors Xmi are formed:(13)Xmi=xi,xi+1,…,xi+m−1−ui,1≤i≤N−m+1,wherem represents embedding dimension and Xmi has m consecutive values. Commencing with the ith point and generalized by removing a baseline:(14)ui=1m∑j=0m−1xi+j.(2)
For eachXmi, the similarity degree between the Xmi and its neighboring vector Xmj is calculated by Di,jm:(15)Di,jm=distdi,jm,r,wheredi,jm denotes the maximum absolute difference in the corresponding scalar components of Xmi and Xmj.(3)
For eachXmi and the fixed tolerance r, let Ai be the number of vectors that satisfy, and then, Bimr is denoted as follows:(16)Bimr=AiN−m+1,1≤i≤N−m.(4)
The average of theBimr is designated as follows:(17)Bmr=1N−m∑i=1N−mBimr,wherer represents the boundary width of the exponential function.(5)
Increasing dimensionm to m+1, steps (1)-(4) are repeated to calculate the corresponding of SE values to find Bm+1r, and SE is defined as follows:(18)SampleEnm,r,N=∑N⟶∞−lnBm+1rBmr.WhenN is finite, the SE can be estimated as follows:(19)SampleEnm,r,N=−lnBm+1rBmr=lnBmr−lnBm+1r.SE determines the time sequence irregularity on the single scale. The smaller the value of FE is, the higher time sequence self-similarity can be achieved. Conversely, the greater the FE value is, the more complicated time sequences without rules can be achieved.


#### 2.2.2. Multiscale Entropy and Multiscale Sample Entropy

ME could be represented as a series of times with different scales, which is obtained through a coarse-grained process. The irregularity of time series is then presented by the ME, where the self-similarity of various scales is reflected. Considering the scale entropy, when there is a sequence with a higher entropy than another sequence, the corresponding complexity will also be higher than the other. That is, if there is a case that the increasing scale of a time series conversely results in the decreased entropy value, there are more chances that the time series holds a relatively simple sequence structure.The MSE is obtained as follows:(1)
Consider a discrete-time seriesXmi:1≤i≤N. Firstly, the value of embedding dimension m and similar tolerance r of the SE are set. Another auxiliary time sequence is then constructed in a vector form, which is presented as ykτ=yk,1τ,yk,2τ,…yk,pτ and named coarse-grained vector.(20)ykτ=1τ∑i=j−1τ+kjτ+k−1Xmi,1≤j≤Nτ,1≤k≤τ,whereτ is the scale factor. Actually, we can find that when τ=1, the coarse-grained time series is equal to the previous one Xmi, which means that by dividing the length of time series, we can decompose the original time series into a coarse-grained vector series ykτ. It is worth mentioning that the number of coarse-grained time series ykτ is τ and the length is N/τ.(2)
Based on the different scalesτ of the time sequence, the SE values are obtained. Generally, the r in SE takes the standard deviation of the primary time series when calculating the value of SE. This procedure is called MFE analysis.(21)MSEX,τ,m,r=SampleEnm,r,y1τ.


#### 2.2.3. Refined Composite Multiscale Sample Entropy

In the process of applying coarse graining, the length of the original time series is declined by introducing a factor ofτ, as shown in equation (13). For example, the authors developed the RCMMSE algorithm with the aim to boost the accuracy of the MSE. In particular, the SE of all the time series is obtained at a certain scale of factor τ. The RCMMSE value is defined as follows:(22)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr,where B¯m+1r=1/τ∑k=1τBm+1r and B¯mr=1/τ∑k=1τBmr. Thus, the RCMMSE value can be obtained as follows:(23)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr=−ln1/τ∑k=1τBm+1r1/τ∑k=1τBmr=−ln∑k=1τBm+1r∑k=1τBmr.


#### 2.2.4. Refined Composite Multivariate Multiscale Sample Entropy

Given time seriesXmpi with p variables, the RCMMSE is calculated according to the following procedures:(1)
Each coarse-grained time seriesykτ=yk,1τ,yk,2τ,…,yk,pτ requires to compute the RCMMSE, the multiple embedding vectors M, and therefore, the original Xmpi is reconstructed using the multiple embedding vectors M.(24)Xmpi=x1i,…x1i+m1−1λ1,x2i,…x2i+m2−1λ2,…x1i,…x1i+m1−1λ1,hereλ=λ1,λ2,…λp is the time delay vector.(2)
The maximum distance of two vectorsXmpi and Xmpi is defined as follows:(25)DXmpi,Xmpj=maxl=1,2,…,mxi+l−1−xj+l−1.(3)
For a composite time delayXmpi and threshold r, the parameter Pi is the number of the vector pairs, and therefore, DXmpi,Xmpj≤r,j≠i, and then, the probability BimrRCMMSE is computed as follows:(26)BimrRCMMSE=PiN−m+1,1≤i≤N−m,herem=maxM∗maxλ.(4)
The mean of theBimrRCMMSE is obtained as follows:(27)B¯imrRCMMSE=1N−m∑i=1N−mBimr.Increasing dimensionm to m+1, steps (1)-(4) are repeated to get the corresponding SE values to find Bim+1rRCMMSE and B¯im+1rRCMMSE.(28)B¯im+1rRCMMSE=1pN−m−1∑i=1pN−m−1Bim+1r.The calculation of the RCMMSE is defined as follows:(29)RCMMSEX,τ,m,r=−lnB¯im+1rRCMMSEB¯imrRCMMSE.


## 2.1. Basic Principle of Variational Mode Decomposition

The VMD process is divided into establishment and solution of variational constraint problems. Assuming there is a limited bandwidth in each mode, the variational problem is formulated ask mode functions ukt, which could minimize the estimated bandwidth. The summation of each mode is constrained to be corresponded to the input signal f. In particular, the structure of this problem could be divided into the following steps.


### 2.1.1. Variational Problem Formulation

①
Hilbert transform: the analytic signal of every mode function is obtained, whose purpose is to get its spectrum:(1)δt+jπt∗ukt.②
The centre frequencye−jωkt of each modal signal is estimated. The spectrum of each mode is modified to be adaptive to its own baseband.(2)δt+jπt∗ukte−jωkt.③
TheL2 is derived, and it denotes the signal gradient, and the signal’s bandwidth is obtained. Therefore, constrained variational questions can be obtained as follows:(3)minuk,ωk∑k∂tδt+jπtμkt−jωkt2s.t.∑kμk=f,where μ is the decomposition of the k mode components, μ=μ1,μ2,…,μk, ω is the centre frequencies of each mode, and ω=ω1,ω2,…,ωk.


### 2.1.2. The Solution of the Problem

①
The above constraint problem could be modified as a nonbinding problem by adding it to the penalty factor∂ and the Lagrangian multiplication operator λt, and here, the second penalty factor is used to guarantee the accuracy of the signal with the Gaussian noise, and the Lagrangian operator is adopted to ensure the constraint condition satisfied strictly, and the extended Lagrangian description is as follows:(4)Luk,ωk,λ≔∂∑k∂tδt+jπt∗utke−jwkt2+ft−∑kukt22+λt,ft−∑kukt.②
VMD used multiplication operator alternating direction method of multipliers (ADMMs) to solve the above variable problem. By alternately updating A, B, and C to seek the “saddle” of the Lagrangian expression:(5)ukn+1=argminuk∈ka∂tδt+jπt∗ukte−jwkt22+ft−∑iuit+λt222,where ωk is equal to ωkn+1, and ∑tuit is equal to ∑i≠kuitn+1.Equation (5) is transformed to another domain by frequency by the Parseval/Plancherel Fourier equidistant transformation:(6)u^kn+1=argminajω1+sgnω+ωk⋅u^kω+ωk22+f^ω−∑iu^iω+λ^ω222.The first item ofω with ω−ωk instead is as follows:(7)u^kn+1=argminajω−ωk1+sgnωu^kω22+f^ω−∑iu^iω+λ^ω222.Equation (7) is converted into a nonnegative integral form with frequency interval:(8)u^kn+1=argminu^k,uk∈X∫0∞4aω−ωk2u^kw2+2f^ω−∑iu^iω+λ^ω22dω.In this case, the solution of the quadratic optimization problem is as follows:(9)u^kn+1ω=f^ω−∑i≠ku^iω+λ^ω/21+2aω−ωk2.Based on the above process, the centre frequency value is converted to the frequency domain:(10)ωkn+1=argmin∫0∞ω−ωk2u^kω2dω.The centre frequency updated scheme is as follows:(11)ωkn+1=∫0∞ωu^kω2dω∫0∞u^kω2dω,where u^kn+1ω is equivalent to the current remaining amount of f^ω−∑i≠ku^iω Wiener filtering. A is the central gravity of the spectrum of the current mode function. For inverse Fourier transmission of u^kω, then the real part is ukt.The procedure of VMD algorithm is as follows:(1)
Initialize theu^k1, ωk1, λ^1, and n.(2)
Updateuk and ωk followed by (9) and (11).(3)
Updateλ, and(12)λ^n+1ω←λ^nω+τf^ω−∑ku^kn+1ω.(4)
For a given accuracye>0, if ∑ku^kn+1−u^kn22/u^kn22<e, stop iteration; otherwise, return to Step 2.


## 2.1.1. Variational Problem Formulation

①
Hilbert transform: the analytic signal of every mode function is obtained, whose purpose is to get its spectrum:(1)δt+jπt∗ukt.②
The centre frequencye−jωkt of each modal signal is estimated. The spectrum of each mode is modified to be adaptive to its own baseband.(2)δt+jπt∗ukte−jωkt.③
TheL2 is derived, and it denotes the signal gradient, and the signal’s bandwidth is obtained. Therefore, constrained variational questions can be obtained as follows:(3)minuk,ωk∑k∂tδt+jπtμkt−jωkt2s.t.∑kμk=f,where μ is the decomposition of the k mode components, μ=μ1,μ2,…,μk, ω is the centre frequencies of each mode, and ω=ω1,ω2,…,ωk.


## 2.1.2. The Solution of the Problem

①
The above constraint problem could be modified as a nonbinding problem by adding it to the penalty factor∂ and the Lagrangian multiplication operator λt, and here, the second penalty factor is used to guarantee the accuracy of the signal with the Gaussian noise, and the Lagrangian operator is adopted to ensure the constraint condition satisfied strictly, and the extended Lagrangian description is as follows:(4)Luk,ωk,λ≔∂∑k∂tδt+jπt∗utke−jwkt2+ft−∑kukt22+λt,ft−∑kukt.②
VMD used multiplication operator alternating direction method of multipliers (ADMMs) to solve the above variable problem. By alternately updating A, B, and C to seek the “saddle” of the Lagrangian expression:(5)ukn+1=argminuk∈ka∂tδt+jπt∗ukte−jwkt22+ft−∑iuit+λt222,where ωk is equal to ωkn+1, and ∑tuit is equal to ∑i≠kuitn+1.Equation (5) is transformed to another domain by frequency by the Parseval/Plancherel Fourier equidistant transformation:(6)u^kn+1=argminajω1+sgnω+ωk⋅u^kω+ωk22+f^ω−∑iu^iω+λ^ω222.The first item ofω with ω−ωk instead is as follows:(7)u^kn+1=argminajω−ωk1+sgnωu^kω22+f^ω−∑iu^iω+λ^ω222.Equation (7) is converted into a nonnegative integral form with frequency interval:(8)u^kn+1=argminu^k,uk∈X∫0∞4aω−ωk2u^kw2+2f^ω−∑iu^iω+λ^ω22dω.In this case, the solution of the quadratic optimization problem is as follows:(9)u^kn+1ω=f^ω−∑i≠ku^iω+λ^ω/21+2aω−ωk2.Based on the above process, the centre frequency value is converted to the frequency domain:(10)ωkn+1=argmin∫0∞ω−ωk2u^kω2dω.The centre frequency updated scheme is as follows:(11)ωkn+1=∫0∞ωu^kω2dω∫0∞u^kω2dω,where u^kn+1ω is equivalent to the current remaining amount of f^ω−∑i≠ku^iω Wiener filtering. A is the central gravity of the spectrum of the current mode function. For inverse Fourier transmission of u^kω, then the real part is ukt.The procedure of VMD algorithm is as follows:(1)
Initialize theu^k1, ωk1, λ^1, and n.(2)
Updateuk and ωk followed by (9) and (11).(3)
Updateλ, and(12)λ^n+1ω←λ^nω+τf^ω−∑ku^kn+1ω.(4)
For a given accuracye>0, if ∑ku^kn+1−u^kn22/u^kn22<e, stop iteration; otherwise, return to Step 2.


## 2.2. Basic Principle of Refined Composite Multivariate Multiscale Sample Entropy

### 2.2.1. Sample Entropy

The SE uses the measurement of the exponential function for two sequences with a tolerancer from m points to remain r of each other at the next point. For a time series with N sample points, the procedure of SE calculation is shown as follows:(1)
Them length vectors Xmi are formed:(13)Xmi=xi,xi+1,…,xi+m−1−ui,1≤i≤N−m+1,wherem represents embedding dimension and Xmi has m consecutive values. Commencing with the ith point and generalized by removing a baseline:(14)ui=1m∑j=0m−1xi+j.(2)
For eachXmi, the similarity degree between the Xmi and its neighboring vector Xmj is calculated by Di,jm:(15)Di,jm=distdi,jm,r,wheredi,jm denotes the maximum absolute difference in the corresponding scalar components of Xmi and Xmj.(3)
For eachXmi and the fixed tolerance r, let Ai be the number of vectors that satisfy, and then, Bimr is denoted as follows:(16)Bimr=AiN−m+1,1≤i≤N−m.(4)
The average of theBimr is designated as follows:(17)Bmr=1N−m∑i=1N−mBimr,wherer represents the boundary width of the exponential function.(5)
Increasing dimensionm to m+1, steps (1)-(4) are repeated to calculate the corresponding of SE values to find Bm+1r, and SE is defined as follows:(18)SampleEnm,r,N=∑N⟶∞−lnBm+1rBmr.WhenN is finite, the SE can be estimated as follows:(19)SampleEnm,r,N=−lnBm+1rBmr=lnBmr−lnBm+1r.SE determines the time sequence irregularity on the single scale. The smaller the value of FE is, the higher time sequence self-similarity can be achieved. Conversely, the greater the FE value is, the more complicated time sequences without rules can be achieved.


### 2.2.2. Multiscale Entropy and Multiscale Sample Entropy

ME could be represented as a series of times with different scales, which is obtained through a coarse-grained process. The irregularity of time series is then presented by the ME, where the self-similarity of various scales is reflected. Considering the scale entropy, when there is a sequence with a higher entropy than another sequence, the corresponding complexity will also be higher than the other. That is, if there is a case that the increasing scale of a time series conversely results in the decreased entropy value, there are more chances that the time series holds a relatively simple sequence structure.The MSE is obtained as follows:(1)
Consider a discrete-time seriesXmi:1≤i≤N. Firstly, the value of embedding dimension m and similar tolerance r of the SE are set. Another auxiliary time sequence is then constructed in a vector form, which is presented as ykτ=yk,1τ,yk,2τ,…yk,pτ and named coarse-grained vector.(20)ykτ=1τ∑i=j−1τ+kjτ+k−1Xmi,1≤j≤Nτ,1≤k≤τ,whereτ is the scale factor. Actually, we can find that when τ=1, the coarse-grained time series is equal to the previous one Xmi, which means that by dividing the length of time series, we can decompose the original time series into a coarse-grained vector series ykτ. It is worth mentioning that the number of coarse-grained time series ykτ is τ and the length is N/τ.(2)
Based on the different scalesτ of the time sequence, the SE values are obtained. Generally, the r in SE takes the standard deviation of the primary time series when calculating the value of SE. This procedure is called MFE analysis.(21)MSEX,τ,m,r=SampleEnm,r,y1τ.


### 2.2.3. Refined Composite Multiscale Sample Entropy

In the process of applying coarse graining, the length of the original time series is declined by introducing a factor ofτ, as shown in equation (13). For example, the authors developed the RCMMSE algorithm with the aim to boost the accuracy of the MSE. In particular, the SE of all the time series is obtained at a certain scale of factor τ. The RCMMSE value is defined as follows:(22)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr,where B¯m+1r=1/τ∑k=1τBm+1r and B¯mr=1/τ∑k=1τBmr. Thus, the RCMMSE value can be obtained as follows:(23)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr=−ln1/τ∑k=1τBm+1r1/τ∑k=1τBmr=−ln∑k=1τBm+1r∑k=1τBmr.


### 2.2.4. Refined Composite Multivariate Multiscale Sample Entropy

Given time seriesXmpi with p variables, the RCMMSE is calculated according to the following procedures:(1)
Each coarse-grained time seriesykτ=yk,1τ,yk,2τ,…,yk,pτ requires to compute the RCMMSE, the multiple embedding vectors M, and therefore, the original Xmpi is reconstructed using the multiple embedding vectors M.(24)Xmpi=x1i,…x1i+m1−1λ1,x2i,…x2i+m2−1λ2,…x1i,…x1i+m1−1λ1,hereλ=λ1,λ2,…λp is the time delay vector.(2)
The maximum distance of two vectorsXmpi and Xmpi is defined as follows:(25)DXmpi,Xmpj=maxl=1,2,…,mxi+l−1−xj+l−1.(3)
For a composite time delayXmpi and threshold r, the parameter Pi is the number of the vector pairs, and therefore, DXmpi,Xmpj≤r,j≠i, and then, the probability BimrRCMMSE is computed as follows:(26)BimrRCMMSE=PiN−m+1,1≤i≤N−m,herem=maxM∗maxλ.(4)
The mean of theBimrRCMMSE is obtained as follows:(27)B¯imrRCMMSE=1N−m∑i=1N−mBimr.Increasing dimensionm to m+1, steps (1)-(4) are repeated to get the corresponding SE values to find Bim+1rRCMMSE and B¯im+1rRCMMSE.(28)B¯im+1rRCMMSE=1pN−m−1∑i=1pN−m−1Bim+1r.The calculation of the RCMMSE is defined as follows:(29)RCMMSEX,τ,m,r=−lnB¯im+1rRCMMSEB¯imrRCMMSE.


## 2.2.1. Sample Entropy

The SE uses the measurement of the exponential function for two sequences with a tolerancer from m points to remain r of each other at the next point. For a time series with N sample points, the procedure of SE calculation is shown as follows:(1)
Them length vectors Xmi are formed:(13)Xmi=xi,xi+1,…,xi+m−1−ui,1≤i≤N−m+1,wherem represents embedding dimension and Xmi has m consecutive values. Commencing with the ith point and generalized by removing a baseline:(14)ui=1m∑j=0m−1xi+j.(2)
For eachXmi, the similarity degree between the Xmi and its neighboring vector Xmj is calculated by Di,jm:(15)Di,jm=distdi,jm,r,wheredi,jm denotes the maximum absolute difference in the corresponding scalar components of Xmi and Xmj.(3)
For eachXmi and the fixed tolerance r, let Ai be the number of vectors that satisfy, and then, Bimr is denoted as follows:(16)Bimr=AiN−m+1,1≤i≤N−m.(4)
The average of theBimr is designated as follows:(17)Bmr=1N−m∑i=1N−mBimr,wherer represents the boundary width of the exponential function.(5)
Increasing dimensionm to m+1, steps (1)-(4) are repeated to calculate the corresponding of SE values to find Bm+1r, and SE is defined as follows:(18)SampleEnm,r,N=∑N⟶∞−lnBm+1rBmr.WhenN is finite, the SE can be estimated as follows:(19)SampleEnm,r,N=−lnBm+1rBmr=lnBmr−lnBm+1r.SE determines the time sequence irregularity on the single scale. The smaller the value of FE is, the higher time sequence self-similarity can be achieved. Conversely, the greater the FE value is, the more complicated time sequences without rules can be achieved.


## 2.2.2. Multiscale Entropy and Multiscale Sample Entropy

ME could be represented as a series of times with different scales, which is obtained through a coarse-grained process. The irregularity of time series is then presented by the ME, where the self-similarity of various scales is reflected. Considering the scale entropy, when there is a sequence with a higher entropy than another sequence, the corresponding complexity will also be higher than the other. That is, if there is a case that the increasing scale of a time series conversely results in the decreased entropy value, there are more chances that the time series holds a relatively simple sequence structure.The MSE is obtained as follows:(1)
Consider a discrete-time seriesXmi:1≤i≤N. Firstly, the value of embedding dimension m and similar tolerance r of the SE are set. Another auxiliary time sequence is then constructed in a vector form, which is presented as ykτ=yk,1τ,yk,2τ,…yk,pτ and named coarse-grained vector.(20)ykτ=1τ∑i=j−1τ+kjτ+k−1Xmi,1≤j≤Nτ,1≤k≤τ,whereτ is the scale factor. Actually, we can find that when τ=1, the coarse-grained time series is equal to the previous one Xmi, which means that by dividing the length of time series, we can decompose the original time series into a coarse-grained vector series ykτ. It is worth mentioning that the number of coarse-grained time series ykτ is τ and the length is N/τ.(2)
Based on the different scalesτ of the time sequence, the SE values are obtained. Generally, the r in SE takes the standard deviation of the primary time series when calculating the value of SE. This procedure is called MFE analysis.(21)MSEX,τ,m,r=SampleEnm,r,y1τ.


## 2.2.3. Refined Composite Multiscale Sample Entropy

In the process of applying coarse graining, the length of the original time series is declined by introducing a factor ofτ, as shown in equation (13). For example, the authors developed the RCMMSE algorithm with the aim to boost the accuracy of the MSE. In particular, the SE of all the time series is obtained at a certain scale of factor τ. The RCMMSE value is defined as follows:(22)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr,where B¯m+1r=1/τ∑k=1τBm+1r and B¯mr=1/τ∑k=1τBmr. Thus, the RCMMSE value can be obtained as follows:(23)RCMMSEX,τ,m,r=−lnB¯m+1rB¯mr=−ln1/τ∑k=1τBm+1r1/τ∑k=1τBmr=−ln∑k=1τBm+1r∑k=1τBmr.


## 2.2.4. Refined Composite Multivariate Multiscale Sample Entropy

Given time seriesXmpi with p variables, the RCMMSE is calculated according to the following procedures:(1)
Each coarse-grained time seriesykτ=yk,1τ,yk,2τ,…,yk,pτ requires to compute the RCMMSE, the multiple embedding vectors M, and therefore, the original Xmpi is reconstructed using the multiple embedding vectors M.(24)Xmpi=x1i,…x1i+m1−1λ1,x2i,…x2i+m2−1λ2,…x1i,…x1i+m1−1λ1,hereλ=λ1,λ2,…λp is the time delay vector.(2)
The maximum distance of two vectorsXmpi and Xmpi is defined as follows:(25)DXmpi,Xmpj=maxl=1,2,…,mxi+l−1−xj+l−1.(3)
For a composite time delayXmpi and threshold r, the parameter Pi is the number of the vector pairs, and therefore, DXmpi,Xmpj≤r,j≠i, and then, the probability BimrRCMMSE is computed as follows:(26)BimrRCMMSE=PiN−m+1,1≤i≤N−m,herem=maxM∗maxλ.(4)
The mean of theBimrRCMMSE is obtained as follows:(27)B¯imrRCMMSE=1N−m∑i=1N−mBimr.Increasing dimensionm to m+1, steps (1)-(4) are repeated to get the corresponding SE values to find Bim+1rRCMMSE and B¯im+1rRCMMSE.(28)B¯im+1rRCMMSE=1pN−m−1∑i=1pN−m−1Bim+1r.The calculation of the RCMMSE is defined as follows:(29)RCMMSEX,τ,m,r=−lnB¯im+1rRCMMSEB¯imrRCMMSE.


## 3. Scheme Comparisons

### 3.1. Comparison of EEMD/VMD

To compare the EEMD and VMD models, a simulation signal is used to prove that the VMD model is better than EEMD. The simulation signalXt is superimposed by the multifrequency signal and the random noise signal with the standard deviation of 1, and the following formula is the signal definition.(30)Xt=1+cos2π×30t×cos2π×125t+1+cos2π×30t×cos2π×155t+cos2π×5t+1+cos2π×30t×cos2π×185t.The time-domain waveforms and its envelope spectrum of the original signal are given in Figure1. As shown in Figure 1(b), the main frequency of the original signal focuses on 0–200 Hz, especially 30.27 Hz. Then, EEMD and VMD were used to decompose the aforementioned signal.Figure 1
Time-domain waveforms and the envelope spectrum of the original signal. (a) Time domain and (b) envelope spectrum.
(a)(b)To achieve the VMD, firstly, the number of modesk should be determined. It should be noted that if the value of k in VMD is chosen too tiny, the original signal with time frequency could not be fully captured by the decomposition of the mode. Larger k values produce the similar frequency between the BL-IMF components, and they may be over-decomposition. Therefore, we use the observing centre frequency of the signal scheme to select the applicable k in this study [19]. The results of the centre frequency corresponding to different k values are given in Table 1.Table 1
Centre frequency under differentk values in VMD.
KCentre frequency (kHz)BL-IMF1BL-IMF2BL-IMF3BL-IMF4BL-IMF5BL-IMF620.00270.0820————30.00220.08830.0621———40.00220.08880.06690.0574——50.00200.08690.06360.05820.0510—60.00200.07980.09430.06030.05100.0513As shown in Table1, there is a centre frequency difference of less than 0.008 kHz from BL-IMF3 to BL-IMF5 modes (0.0636 kHz in BL-IMF3, 0.0582 kHz in BL-IMF4, and 0.051 kHz in BL-IMF5) when k=5. These three modes exist similar frequency. Therefore, the number of k values is selected as 4 in this study to decompose the above simulation signal.The EEMD/VMD models were used to decompose the original signal, and the results of IMFs (EEMD) and BL-IMFs (VMD) are given in Figures1(a) and 2(a). Figures 1(b) and 2(b) are the corresponding spectrum analysis results. As shown in Figure 2, each BL-IMF component is mainly distributed around a single frequency (30.27 Hz), compared with VMD in Figure 1(b), and some IMF components in Figure 2(b) have a series frequency from IMF4 mode.Figure 2
VMD signal decomposition. (a) The BL-IMF components. (b) The spectrum of each BL-IMF.
(a)(b)The mentioned result shows that the decomposition effect of EEMD algorithm is not ideal for the multicomponent synthesis simulation signal, and the mode mixing is serious. Because some slight signals drown in the signal, which is to be decomposed, EEMD in the selection process of the three-spline envelope fitting leads to decomposition bias. The weak signal is embedded in the vast majority of the strong signal where the EEMD can be filtered and extracted, but when the weak signal appears only in the maximum slope range of the strong signal, the weak signal will be in the form of wave frequency modulation and does not produce additional local extreme points. Therefore, EEMD is difficult to extract the useful components and easy to produce some components with mode mixing.It can be seen from Figure1(b) that VMD can not only effectively remove the dummy components, but also each BL-IMF exhibits a mode in the range of a certain scale, and there is no mode mixing problem between them. Therefore, the scale characterization and decomposition effect of VMD and better than EEMD, and this indicates that VMD has good robustness Figure 3.Figure 3
Signal decomposition (EEMD). (a) The IMF components. (b) The spectrum of each IMF mode.
(a)(b)


### 3.2. Comparison of Validity of MSE/RCMMSE

In this section, the1/f noise signal was employed to demonstrate the superior performance of the RCMMSE compared with MSE. Resorting to MSE and RCMMSE, the probabilities of inducing undefined entropy were presented to address the white noise and 1/f noises with 200 samples. For each sample, the length is chosen as 1000. The results of the corresponding probabilities are shown in Table 2. From Table 2, we can see that the probability of undefined entropy showed the same tendency as the time scale increased, while there was a decrease in the length of the time series. As given in Table 1, the probability of undefined entropy is about zero when they were used to analyse white noise in the MSE model. On the contrary, the probability of inducing undefined entropy is 0.01 when τ=4 for all 200 1/f noise samples. In the MSE algorithm, the entropy is undefined when Bmr or Bm+1r is zero. However, the RCMMSE can successfully get the values of entropy from 1 to 20 when the white and 1/f noises are considered. Hence, the RCMMSE model is superior to MSE.Table 2
Probabilities of undefined entropy.
Scale factorMSERCMMSEWhite noise1/f noiseWhite noise1/f noise1–30000400.0100500.01500600.04500700.1100800.15500900.275001000.35001100.36001200.4824001300.4949001400.5606001500.5377001600.5714001700.7059001800.7211001900.6402002000.666700


### 3.3. Comparison of Accuracy in MSE/RCMMSE

In this section, the efficiency of the MSE/RCMMSE models was verified through case studies. In each simulation study, we employ noise samples with 1000 data points containing white and1/f noises. The results of MSE/RCMMSEs are presented in Figure 4. Meanwhile, Figures 5 and 6 present the results of means and standard deviations of entropies.Figure 4
MER/RCMMSE values of white and 1/f noises. (a) White noise. (b) 1/f noise.
(a)(b)Figure 5
Mean deviation of 200 analysis results of white and 1/f noises using MSE and RCMMSE. (a) White noise. (b) 1/f noise.
(a)(b)Figure 6
Standard deviation of 200 analysis results of white noise and 1/f noises using MSE and RCMMSE. (a) White noise. (b) 1/f noise.
(a)(b)As seen in Figure4(a), the overall recently RCMMSE model shows higher performance compared with the MSE model in the condition of white noise, resulting from the tolerance r in SE, which is utilized to evaluate the similarity between any two time series. Note that r was often chosen as 0.15 × SD of the original time series.It can be seen from Figures5(a) and 5(b) that the means of the entropy values obtained using the MSE/RCMMSE are nearly equal in white noise. Nevertheless, the means of the entropy values of the MSE are higher than that of the RCMMSE (see Figure 5(b)). Figure 6 shows that the standard deviation of RCMMSE is all lower than that of MSE, and this result indicates that the entropies obtained using the RCMMSE algorithm were more consistent than those obtained using the MSE algorithm.With the purpose of studying the relationship between the data length and the effectiveness of the MSE/RCMMSE models, several different data lengths (N = 600, 1200, 1800, and 2400) are used to get the entropies. In Table 3, the means and standard deviation of the entropies are presented. It can be found that the results of the means of entropies by MSE and RCMMSE are nearly equivalent with the long primary time series (N = 2400). However, with short original time series, the results of means show a significant difference. On the other hand, the results of the RCMMSE algorithm hold the lowest error when the discrimination between the entropy of 1/f noise with 500 data points and the entropy of 1/f noise with 2000 data points is considered. That is to say, the RCMMSE has higher performance compared with other algorithms at any length of the time series. In conclusion, the RCMMSE is more reliable than MSE according to the above discussions.Table 3
Means and standard deviations of the MSE and RCMMSE of 1/f noise at a scale factor of 20.
Data lengthMSERCMMSEMeansStandard deviationsMeansStandard deviations500——2.13230.52341000——2.01850.26671500——2.01850.17612000——1.93730.1230


## 3.1. Comparison of EEMD/VMD

To compare the EEMD and VMD models, a simulation signal is used to prove that the VMD model is better than EEMD. The simulation signalXt is superimposed by the multifrequency signal and the random noise signal with the standard deviation of 1, and the following formula is the signal definition.(30)Xt=1+cos2π×30t×cos2π×125t+1+cos2π×30t×cos2π×155t+cos2π×5t+1+cos2π×30t×cos2π×185t.The time-domain waveforms and its envelope spectrum of the original signal are given in Figure1. As shown in Figure 1(b), the main frequency of the original signal focuses on 0–200 Hz, especially 30.27 Hz. Then, EEMD and VMD were used to decompose the aforementioned signal.Figure 1
Time-domain waveforms and the envelope spectrum of the original signal. (a) Time domain and (b) envelope spectrum.
(a)(b)To achieve the VMD, firstly, the number of modesk should be determined. It should be noted that if the value of k in VMD is chosen too tiny, the original signal with time frequency could not be fully captured by the decomposition of the mode. Larger k values produce the similar frequency between the BL-IMF components, and they may be over-decomposition. Therefore, we use the observing centre frequency of the signal scheme to select the applicable k in this study [19]. The results of the centre frequency corresponding to different k values are given in Table 1.Table 1
Centre frequency under differentk values in VMD.
KCentre frequency (kHz)BL-IMF1BL-IMF2BL-IMF3BL-IMF4BL-IMF5BL-IMF620.00270.0820————30.00220.08830.0621———40.00220.08880.06690.0574——50.00200.08690.06360.05820.0510—60.00200.07980.09430.06030.05100.0513As shown in Table1, there is a centre frequency difference of less than 0.008 kHz from BL-IMF3 to BL-IMF5 modes (0.0636 kHz in BL-IMF3, 0.0582 kHz in BL-IMF4, and 0.051 kHz in BL-IMF5) when k=5. These three modes exist similar frequency. Therefore, the number of k values is selected as 4 in this study to decompose the above simulation signal.The EEMD/VMD models were used to decompose the original signal, and the results of IMFs (EEMD) and BL-IMFs (VMD) are given in Figures1(a) and 2(a). Figures 1(b) and 2(b) are the corresponding spectrum analysis results. As shown in Figure 2, each BL-IMF component is mainly distributed around a single frequency (30.27 Hz), compared with VMD in Figure 1(b), and some IMF components in Figure 2(b) have a series frequency from IMF4 mode.Figure 2
VMD signal decomposition. (a) The BL-IMF components. (b) The spectrum of each BL-IMF.
(a)(b)The mentioned result shows that the decomposition effect of EEMD algorithm is not ideal for the multicomponent synthesis simulation signal, and the mode mixing is serious. Because some slight signals drown in the signal, which is to be decomposed, EEMD in the selection process of the three-spline envelope fitting leads to decomposition bias. The weak signal is embedded in the vast majority of the strong signal where the EEMD can be filtered and extracted, but when the weak signal appears only in the maximum slope range of the strong signal, the weak signal will be in the form of wave frequency modulation and does not produce additional local extreme points. Therefore, EEMD is difficult to extract the useful components and easy to produce some components with mode mixing.It can be seen from Figure1(b) that VMD can not only effectively remove the dummy components, but also each BL-IMF exhibits a mode in the range of a certain scale, and there is no mode mixing problem between them. Therefore, the scale characterization and decomposition effect of VMD and better than EEMD, and this indicates that VMD has good robustness Figure 3.Figure 3
Signal decomposition (EEMD). (a) The IMF components. (b) The spectrum of each IMF mode.
(a)(b)


## 3.2. Comparison of Validity of MSE/RCMMSE

In this section, the1/f noise signal was employed to demonstrate the superior performance of the RCMMSE compared with MSE. Resorting to MSE and RCMMSE, the probabilities of inducing undefined entropy were presented to address the white noise and 1/f noises with 200 samples. For each sample, the length is chosen as 1000. The results of the corresponding probabilities are shown in Table 2. From Table 2, we can see that the probability of undefined entropy showed the same tendency as the time scale increased, while there was a decrease in the length of the time series. As given in Table 1, the probability of undefined entropy is about zero when they were used to analyse white noise in the MSE model. On the contrary, the probability of inducing undefined entropy is 0.01 when τ=4 for all 200 1/f noise samples. In the MSE algorithm, the entropy is undefined when Bmr or Bm+1r is zero. However, the RCMMSE can successfully get the values of entropy from 1 to 20 when the white and 1/f noises are considered. Hence, the RCMMSE model is superior to MSE.Table 2
Probabilities of undefined entropy.
Scale factorMSERCMMSEWhite noise1/f noiseWhite noise1/f noise1–30000400.0100500.01500600.04500700.1100800.15500900.275001000.35001100.36001200.4824001300.4949001400.5606001500.5377001600.5714001700.7059001800.7211001900.6402002000.666700


## 3.3. Comparison of Accuracy in MSE/RCMMSE

In this section, the efficiency of the MSE/RCMMSE models was verified through case studies. In each simulation study, we employ noise samples with 1000 data points containing white and1/f noises. The results of MSE/RCMMSEs are presented in Figure 4. Meanwhile, Figures 5 and 6 present the results of means and standard deviations of entropies.Figure 4
MER/RCMMSE values of white and 1/f noises. (a) White noise. (b) 1/f noise.
(a)(b)Figure 5
Mean deviation of 200 analysis results of white and 1/f noises using MSE and RCMMSE. (a) White noise. (b) 1/f noise.
(a)(b)Figure 6
Standard deviation of 200 analysis results of white noise and 1/f noises using MSE and RCMMSE. (a) White noise. (b) 1/f noise.
(a)(b)As seen in Figure4(a), the overall recently RCMMSE model shows higher performance compared with the MSE model in the condition of white noise, resulting from the tolerance r in SE, which is utilized to evaluate the similarity between any two time series. Note that r was often chosen as 0.15 × SD of the original time series.It can be seen from Figures5(a) and 5(b) that the means of the entropy values obtained using the MSE/RCMMSE are nearly equal in white noise. Nevertheless, the means of the entropy values of the MSE are higher than that of the RCMMSE (see Figure 5(b)). Figure 6 shows that the standard deviation of RCMMSE is all lower than that of MSE, and this result indicates that the entropies obtained using the RCMMSE algorithm were more consistent than those obtained using the MSE algorithm.With the purpose of studying the relationship between the data length and the effectiveness of the MSE/RCMMSE models, several different data lengths (N = 600, 1200, 1800, and 2400) are used to get the entropies. In Table 3, the means and standard deviation of the entropies are presented. It can be found that the results of the means of entropies by MSE and RCMMSE are nearly equivalent with the long primary time series (N = 2400). However, with short original time series, the results of means show a significant difference. On the other hand, the results of the RCMMSE algorithm hold the lowest error when the discrimination between the entropy of 1/f noise with 500 data points and the entropy of 1/f noise with 2000 data points is considered. That is to say, the RCMMSE has higher performance compared with other algorithms at any length of the time series. In conclusion, the RCMMSE is more reliable than MSE according to the above discussions.Table 3
Means and standard deviations of the MSE and RCMMSE of 1/f noise at a scale factor of 20.
Data lengthMSERCMMSEMeansStandard deviationsMeansStandard deviations500——2.13230.52341000——2.01850.26671500——2.01850.17612000——1.93730.1230


## 4. Experimental Setup

### 4.1. The Dataset Source

The performance of the proposed scheme is verified through experiments in this part. Case Western Reserve University Bearing Data are used, which was obtained from an induction motor. It should be noted that the data are generated by accelerometers at the drive end and fan end. Besides, EDM is employed to model the possible faults in the motor bearings. In particular, the fault tolerance is set as 0.1778 mm. After the fault detection, the data located in the outer race were aggregated in a state, which holds the time configuration at 6:00clock. To ensure the data acquisition system suitable for the vibration signals, the amplifier is particularly designed with a high bandwidth. On the other hand, to improve the accuracy, a sampling frequency of 12000 Hz of each channel was set for the data recorder.Table4 demonstrates different working conditions, which are considered in this experiment. Note that 0.1778 mm is the fault tolerance. “Normal” and “slight” denote the fault severity. In addition, considering the drive end of the motor 1785 rpm is chosen as the motor revolving speed. In total, 200 data samples are adopted with each fault condition consisting of 51 samples. In particular, in each data sample, there are 2048 data points.Table 4
Experimental data of the roller bearings.
Fault categoryFault tolerance (mm)Motor speed (rpm)Number of samplesThe fault severityNR0178551NormalIRF0.1778178551SlightBF0.1778178551SlightORF0.1778178551Slight


### 4.2. Procedure of Our Proposed Scheme

The procedure of the proposed scheme can be shown as follows:(1)
Resorting to EEMD and VMD models, the vibration signals under different cases were disintegrated into a sequence of IMF and BL-IMF modes.(2)
Using the RCMMSE model to get the entropy values, the entropy values were regarded as the input of SVM/RF/PNN models for training and testing.(3)
Testing samples were regarded as the input of the trained SVM/RF/PNN classifier, while the operating conditions can be recognized by the output of the SVM/RF/PNN classifiers. The procedure of our proposed fault diagnosis scheme is shown in Figure7.Figure 7
Procedure of our proposed scheme.


### 4.3. Parameter Selection

There are some parameters in different models that should be present before its calculation:(1)
EEMD: considering EEMD, there are two parameters that need to be determined, namely the ensemble numberm and nit related to the white noise. In particular, nit denotes the amplitude. It should be mentioned that there is a one-to-one correspondence between the result and the ensemble number. Note that the ensemble number should be chosen as a few hundred. In addition, consider the situation that the input signal holds the standard deviation, and a fraction of one percent of error will be caused by the remaining noise with the added noise. In particular, m is set as 100.(2)
VMD: the number of BL-IMF modesk is according to the centre frequency. The second penalty factor ∂ in equation (4) is often set as 2000.(3)
SE/RCMMSE:①
SE: before the process of SE, there are three parameters that must be chosen suitably. Firstly, considering the different lengths of sequences in SE, embedding dimensionm should be designed. In particular, the dynamic process will be reconstructed in a more detailed manner with a larger m. However, it should also be mentioned that if m is chosen too large, the need of N=10m−30m would be greatly limited, which will cause the loss of important information, and the general condition will be hard to satisfy. In the literature, m is often set as 2. Another two parameters, namely, similarity tolerance n and nr, are dependent on the gradient of the exponential function and its bound, respectively. Experimentally, r=0.1−0.25SD, and r is chosen as 0.15SD in SE.②
RCMMSE: the multiple embedding vector in RCMMSE scheme is set asM=21,22,23,…,2k, and here, k is the number of the BL-IMF components, and therefore, the time delay is set as 1.(4)
SVM: the kernel function in SVM is selected as the radial basis function (RBF).(5)
PNN: the distribution density of PNN is set as 1.5.(6)
RF: there are two parameters that need to be confirmed before using RF model, such as the number of input variablesmtry is selected randomly based on the M input variables. Regarding the scale factor τ as 20 and EEMD/VMD-RCMMSE as feature, the number of input variables M=τ=20, and the parameters often meet the condition mtry≤M [15]. Therefore, we let mtry=4 and the number of the DT is set as 4 and 1000 in this simulation.


## 4.1. The Dataset Source

The performance of the proposed scheme is verified through experiments in this part. Case Western Reserve University Bearing Data are used, which was obtained from an induction motor. It should be noted that the data are generated by accelerometers at the drive end and fan end. Besides, EDM is employed to model the possible faults in the motor bearings. In particular, the fault tolerance is set as 0.1778 mm. After the fault detection, the data located in the outer race were aggregated in a state, which holds the time configuration at 6:00clock. To ensure the data acquisition system suitable for the vibration signals, the amplifier is particularly designed with a high bandwidth. On the other hand, to improve the accuracy, a sampling frequency of 12000 Hz of each channel was set for the data recorder.Table4 demonstrates different working conditions, which are considered in this experiment. Note that 0.1778 mm is the fault tolerance. “Normal” and “slight” denote the fault severity. In addition, considering the drive end of the motor 1785 rpm is chosen as the motor revolving speed. In total, 200 data samples are adopted with each fault condition consisting of 51 samples. In particular, in each data sample, there are 2048 data points.Table 4
Experimental data of the roller bearings.
Fault categoryFault tolerance (mm)Motor speed (rpm)Number of samplesThe fault severityNR0178551NormalIRF0.1778178551SlightBF0.1778178551SlightORF0.1778178551Slight


## 4.2. Procedure of Our Proposed Scheme

The procedure of the proposed scheme can be shown as follows:(1)
Resorting to EEMD and VMD models, the vibration signals under different cases were disintegrated into a sequence of IMF and BL-IMF modes.(2)
Using the RCMMSE model to get the entropy values, the entropy values were regarded as the input of SVM/RF/PNN models for training and testing.(3)
Testing samples were regarded as the input of the trained SVM/RF/PNN classifier, while the operating conditions can be recognized by the output of the SVM/RF/PNN classifiers. The procedure of our proposed fault diagnosis scheme is shown in Figure7.Figure 7
Procedure of our proposed scheme.


## 4.3. Parameter Selection

There are some parameters in different models that should be present before its calculation:(1)
EEMD: considering EEMD, there are two parameters that need to be determined, namely the ensemble numberm and nit related to the white noise. In particular, nit denotes the amplitude. It should be mentioned that there is a one-to-one correspondence between the result and the ensemble number. Note that the ensemble number should be chosen as a few hundred. In addition, consider the situation that the input signal holds the standard deviation, and a fraction of one percent of error will be caused by the remaining noise with the added noise. In particular, m is set as 100.(2)
VMD: the number of BL-IMF modesk is according to the centre frequency. The second penalty factor ∂ in equation (4) is often set as 2000.(3)
SE/RCMMSE:①
SE: before the process of SE, there are three parameters that must be chosen suitably. Firstly, considering the different lengths of sequences in SE, embedding dimensionm should be designed. In particular, the dynamic process will be reconstructed in a more detailed manner with a larger m. However, it should also be mentioned that if m is chosen too large, the need of N=10m−30m would be greatly limited, which will cause the loss of important information, and the general condition will be hard to satisfy. In the literature, m is often set as 2. Another two parameters, namely, similarity tolerance n and nr, are dependent on the gradient of the exponential function and its bound, respectively. Experimentally, r=0.1−0.25SD, and r is chosen as 0.15SD in SE.②
RCMMSE: the multiple embedding vector in RCMMSE scheme is set asM=21,22,23,…,2k, and here, k is the number of the BL-IMF components, and therefore, the time delay is set as 1.(4)
SVM: the kernel function in SVM is selected as the radial basis function (RBF).(5)
PNN: the distribution density of PNN is set as 1.5.(6)
RF: there are two parameters that need to be confirmed before using RF model, such as the number of input variablesmtry is selected randomly based on the M input variables. Regarding the scale factor τ as 20 and EEMD/VMD-RCMMSE as feature, the number of input variables M=τ=20, and the parameters often meet the condition mtry≤M [15]. Therefore, we let mtry=4 and the number of the DT is set as 4 and 1000 in this simulation.


## 5. Experimental Results and Analysis

### 5.1. EEMD/VMD

The data are selected from the experiments in which SKF bearings are used. The number of each sample is set as NR:1–50, IRF: 51–100, BF:101–150, ad ORF:151–200. As limited space, here with a sample of each state for an example, the time-domain waveforms of vibration signals under different working conditions are shown in Figure8.Figure 8
Time-domain waveforms of each working condition.The vertical axis is the acceleration vibration amplitude. Because of the influence of noise, it is difficult to find significant differences in different states. As shown in Figure8, it is hard to distinguish the four signals; in particular, there is no obvious regularity in two states of NR and BF signals.The EEMD and VMD were used to disintegrate the vibration signals into series modes (IMF and BL-IMF). Here, we use the IRF vibration signals to observe the frequency change situation, and the waveforms of one IRF vibration signal and its spectrum envelope are shown in Figure9. As shown in Figure 5(b), the IRF signal fault frequency is 164.1 Hz. The results of the centre frequency under different k values with VMD are given in Table 5. There are some modes, such as BL-IMF4 (0.2787 kHz) and BL-IMF5 (0.3071 kHz), and with similar centre frequencies, when k=5 in Table 5, it is considered that there is over-decomposition. Therefore, the number of k values is selected as 4 to decompose the vibration signals.Figure 9
Waveforms and the envelope spectrum of the original signal in time domain. (a) Time domain. (b) Envelope spectrum.
(a)(b)Table 5
Centre frequency under differentk in VMD.
KCentral frequency (KHz)BL-IMF1BL-IMF2BL-IMF3BL-IMF4BL-IMF5BL-IMF620.09220.2631————30.08950.22460.2873———40.04940.11400.22680.2878——50.04900.11370.22390.27870.3071—60.04890.11360.22280.27300.28900.3144Employing EEMD/VMD, the roller bearing primary signals in Table4 are decomposed into IMF1-10 and BL-IMF1-4. Figures 10 and 11 show the VMD and the corresponding spectrum of each mode for an IRF signal.Figure 10
Signal decomposition based on VMD. (a) The BL-IMF components. (b) The spectrum of every BL-IMF.
(a)(b)Figure 11
NR signal based on EEMD. (a) The IMF components. (b) The spectrum of each IMF.
(a)(b)It can be shown from Figure6(b) that the VMD can not only effectively remove the dummy components, but also each IMF in the range of a certain scale, and there is no mode mixing problem between them. Therefore, the scale characterization and decomposition effect of VMD are better than EEMD, and this result reveals that the VMD owns a good robustness.The correlation coefficient scheme is utilized to prove the degree of correlation between each mode and the primary signal to further prove the effectiveness of VMD scheme. In Tables6 and 7, the overall average correlation values of BL-IMF component and original IRF signals are higher than IMF. Therefore, the decomposition effect of VMD is superior to EEMD.Table 6
Total average correlation values of each IMF with EEMD.
ModeIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10EEMDNR0.540.620.460.450.450.270.120.0220.00470.0089IRF0.910.380.1780.1010.0520.02010.001490.000860.000252.46e−05BF0.840.4060.3280.2660.0850.0210.0040.0030.00130.0013ORF0.860.3080.3660.1890.1250.0530.00410.002680.00130.0001Table 7
Total average correlation values of each BL-IMF with VMD.
ModeModeBL-IMF1BL-IMF2BL-IMF3BL-IMF4VMD (k = 4)NR0.61520.44630.49400.5384IRF0.22750.37180.57730.7702BF0.38490.41410.49430.7376ORF0.42010.51870.74210.2708After the vibration signal feature extraction using VMD and EEMD, the RCMMSE was employed to get the entropy value. The value of VMD-RCMMSE values under differentk values and scale factors is shown in Figure 12.Figure 12
RCMMSE values under differentk values and scale factors.Figure12 shows that the RCMMSE values of the VMD effective mode decrease with the increasement in the scale factor and the number of modes k, and the overall RCMMSE values are close to the steady state when the value of k is more than 3, because the VMD algorithm decomposes the original signals into several BL-IMFs with information of different frequency bands. This can also be supported from Table 5. The correlation coefficient of most EEMD modes and the original signal is much smaller than that of VMD mode. This indicates that the VMD algorithm can overcome the mode mixing shortcoming in EEMD. The entropy values of EEMD-RCMMSE and VMD-RCMMSE (k = 4) are shown in Figure 13.Figure 13
RCMMSE values with EEMD and VMD (k = 4) modes.
(a)(b)As shown in Figure13(a), it is difficult to distinguish the four types of roller bearing signals with EEMD-RCMMSE, especially the NR and IRF, BF, and ORF. On the one hand, the complexity of signals under normal states and the conditions with ball fault is similar. In other words, the feature of the two kinds of signals is alike. On the other hand, much noise hidden in the signals may affect the identification of them. However, the RCMMSE curves of BL-IMF components can recognize the four states clearly in Figure 13(b).


### 5.2. Fault Identification

Regarding the extracted RCMMSE vectors as the inputs of the chosen fault classifier, such as SVM/RF/PNN, there are four kinds of data with 200 samples in total. 10, 20, and 30 samples are used under different statuses as the training samples, and the rest samples (20, 30, 40) as the testing samples are used to compare the classification accuracy. Part of the experiment classification results is shown in Figure14 and Table 8.Figure 14
Waveforms in the time domain.
(a)(b)(c)(d)(e)(f)Table 8
Mean and standard deviation of the MSE/MFE/RCMSE and RCMMSE of 1/f noise at a scale factor of 20.
ModePrecision (%)Total number of testing samples16012080EEMD-RCMMSE-SVM69.37573.3386.25EEMD-RCMMSE-RF63.7569.1780EEMD-RCMMSE-PNN71.2574.1671.25VMD-RCMMSE-SVM (k = 3)8588.391.25VMD-RCMMSE-RF (k = 3)81.2588.391.25VMD-RCMMSE-PNN (k = 3)83.1258586.25VMD-RCMMSE-SVM (k = 4)98.12599.1667100VMD-RCMMSE-RF (k = 4)96.2596.6797.5VMD-RCMMSE-PNN (k = 4)93.7590.8391.25VMD-RCMMSE-SVM (k = 5)98.7597.598.75VMD-RCMMSE-RF (k = 5)97.59596.25VMD-RCMMSE-PNN (k = 5)86.8759091.25As given in Table8, the highest classification precision is up to 100% when k=4 in VMD-RCMMSE model. On the contrary, the lowest accuracy is 63.75% in EEMD-RCMMSE model. The overall classification accuracy in VMD-RCMMSE model is higher than EEMD-RCMMSE model. The overall classification accuracy in SVM model is higher than RF/PNN models. As mentioned above, the experimental results show that VMD-RCMMSE combination model can recognize different working conditions of the roller bearings effectively.


## 5.1. EEMD/VMD

The data are selected from the experiments in which SKF bearings are used. The number of each sample is set as NR:1–50, IRF: 51–100, BF:101–150, ad ORF:151–200. As limited space, here with a sample of each state for an example, the time-domain waveforms of vibration signals under different working conditions are shown in Figure8.Figure 8
Time-domain waveforms of each working condition.The vertical axis is the acceleration vibration amplitude. Because of the influence of noise, it is difficult to find significant differences in different states. As shown in Figure8, it is hard to distinguish the four signals; in particular, there is no obvious regularity in two states of NR and BF signals.The EEMD and VMD were used to disintegrate the vibration signals into series modes (IMF and BL-IMF). Here, we use the IRF vibration signals to observe the frequency change situation, and the waveforms of one IRF vibration signal and its spectrum envelope are shown in Figure9. As shown in Figure 5(b), the IRF signal fault frequency is 164.1 Hz. The results of the centre frequency under different k values with VMD are given in Table 5. There are some modes, such as BL-IMF4 (0.2787 kHz) and BL-IMF5 (0.3071 kHz), and with similar centre frequencies, when k=5 in Table 5, it is considered that there is over-decomposition. Therefore, the number of k values is selected as 4 to decompose the vibration signals.Figure 9
Waveforms and the envelope spectrum of the original signal in time domain. (a) Time domain. (b) Envelope spectrum.
(a)(b)Table 5
Centre frequency under differentk in VMD.
KCentral frequency (KHz)BL-IMF1BL-IMF2BL-IMF3BL-IMF4BL-IMF5BL-IMF620.09220.2631————30.08950.22460.2873———40.04940.11400.22680.2878——50.04900.11370.22390.27870.3071—60.04890.11360.22280.27300.28900.3144Employing EEMD/VMD, the roller bearing primary signals in Table4 are decomposed into IMF1-10 and BL-IMF1-4. Figures 10 and 11 show the VMD and the corresponding spectrum of each mode for an IRF signal.Figure 10
Signal decomposition based on VMD. (a) The BL-IMF components. (b) The spectrum of every BL-IMF.
(a)(b)Figure 11
NR signal based on EEMD. (a) The IMF components. (b) The spectrum of each IMF.
(a)(b)It can be shown from Figure6(b) that the VMD can not only effectively remove the dummy components, but also each IMF in the range of a certain scale, and there is no mode mixing problem between them. Therefore, the scale characterization and decomposition effect of VMD are better than EEMD, and this result reveals that the VMD owns a good robustness.The correlation coefficient scheme is utilized to prove the degree of correlation between each mode and the primary signal to further prove the effectiveness of VMD scheme. In Tables6 and 7, the overall average correlation values of BL-IMF component and original IRF signals are higher than IMF. Therefore, the decomposition effect of VMD is superior to EEMD.Table 6
Total average correlation values of each IMF with EEMD.
ModeIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10EEMDNR0.540.620.460.450.450.270.120.0220.00470.0089IRF0.910.380.1780.1010.0520.02010.001490.000860.000252.46e−05BF0.840.4060.3280.2660.0850.0210.0040.0030.00130.0013ORF0.860.3080.3660.1890.1250.0530.00410.002680.00130.0001Table 7
Total average correlation values of each BL-IMF with VMD.
ModeModeBL-IMF1BL-IMF2BL-IMF3BL-IMF4VMD (k = 4)NR0.61520.44630.49400.5384IRF0.22750.37180.57730.7702BF0.38490.41410.49430.7376ORF0.42010.51870.74210.2708After the vibration signal feature extraction using VMD and EEMD, the RCMMSE was employed to get the entropy value. The value of VMD-RCMMSE values under differentk values and scale factors is shown in Figure 12.Figure 12
RCMMSE values under differentk values and scale factors.Figure12 shows that the RCMMSE values of the VMD effective mode decrease with the increasement in the scale factor and the number of modes k, and the overall RCMMSE values are close to the steady state when the value of k is more than 3, because the VMD algorithm decomposes the original signals into several BL-IMFs with information of different frequency bands. This can also be supported from Table 5. The correlation coefficient of most EEMD modes and the original signal is much smaller than that of VMD mode. This indicates that the VMD algorithm can overcome the mode mixing shortcoming in EEMD. The entropy values of EEMD-RCMMSE and VMD-RCMMSE (k = 4) are shown in Figure 13.Figure 13
RCMMSE values with EEMD and VMD (k = 4) modes.
(a)(b)As shown in Figure13(a), it is difficult to distinguish the four types of roller bearing signals with EEMD-RCMMSE, especially the NR and IRF, BF, and ORF. On the one hand, the complexity of signals under normal states and the conditions with ball fault is similar. In other words, the feature of the two kinds of signals is alike. On the other hand, much noise hidden in the signals may affect the identification of them. However, the RCMMSE curves of BL-IMF components can recognize the four states clearly in Figure 13(b).


## 5.2. Fault Identification

Regarding the extracted RCMMSE vectors as the inputs of the chosen fault classifier, such as SVM/RF/PNN, there are four kinds of data with 200 samples in total. 10, 20, and 30 samples are used under different statuses as the training samples, and the rest samples (20, 30, 40) as the testing samples are used to compare the classification accuracy. Part of the experiment classification results is shown in Figure14 and Table 8.Figure 14
Waveforms in the time domain.
(a)(b)(c)(d)(e)(f)Table 8
Mean and standard deviation of the MSE/MFE/RCMSE and RCMMSE of 1/f noise at a scale factor of 20.
ModePrecision (%)Total number of testing samples16012080EEMD-RCMMSE-SVM69.37573.3386.25EEMD-RCMMSE-RF63.7569.1780EEMD-RCMMSE-PNN71.2574.1671.25VMD-RCMMSE-SVM (k = 3)8588.391.25VMD-RCMMSE-RF (k = 3)81.2588.391.25VMD-RCMMSE-PNN (k = 3)83.1258586.25VMD-RCMMSE-SVM (k = 4)98.12599.1667100VMD-RCMMSE-RF (k = 4)96.2596.6797.5VMD-RCMMSE-PNN (k = 4)93.7590.8391.25VMD-RCMMSE-SVM (k = 5)98.7597.598.75VMD-RCMMSE-RF (k = 5)97.59596.25VMD-RCMMSE-PNN (k = 5)86.8759091.25As given in Table8, the highest classification precision is up to 100% when k=4 in VMD-RCMMSE model. On the contrary, the lowest accuracy is 63.75% in EEMD-RCMMSE model. The overall classification accuracy in VMD-RCMMSE model is higher than EEMD-RCMMSE model. The overall classification accuracy in SVM model is higher than RF/PNN models. As mentioned above, the experimental results show that VMD-RCMMSE combination model can recognize different working conditions of the roller bearings effectively.


## 6. Conclusion

The VMD analysis model and RCMMSE are proposed in this study. The VMD model is utilized to analyse the complex vibration signals. Therefore, BL-IMF can be obtained from the decomposition of the roller bearing vibration signals, the RCMMSE can deal with multichannel data, and the EEMD/VMD-RCMMMSE values are regarded as eigenvectors. Because the VMD can achieve adaptive subdivision of each component in the frequency domain of signals, the roller bearing vibration signals are distinguished. Finally, the VMD-RCMMSE-SVM/RF/PNN combination models can distinguish roller bearing fault types effectively.With the rapid development of Internet of things (IoT) [1, 2] and Industry 4.0 [3], there are increasingly massive real-time data from various types of mechanical equipment [4]. The availability of these data that contain abundant information about machine health has attracted more and more enterprises’ attention. It has been proved that large volume, high velocity, and diversity mechanical big data are the major properties of mechanical big data [5, 6]. Effective feature extraction from these data and accurate machinery health state evaluation with ever-accelerated updating of schemes have become hot research issues in the prognostic and health management systems in the era of industrial IoT [7].


## 7. Future Ideas

Most related roller bearing signal fault diagnoses depend on manual aids and expert knowledge. To diagnose more smartly, convolutional neural network (\textsc{CNN}), as well as long short-term memory (LSTM), will be considered for an end-to-end scenario. For ensuring high recognition accuracy, the fault classifier based on CNN may be used to extract the effective signal from the severely distorted signal; the LSTM could cope with time-varying bearing failure. Additionally, the aforementioned SVM fault classifier performs well when the data samples are insufficient. In the future, we may consider the generative adversarial networks (GANs) and transfer learning for lack of data. By simulating the distribution of fault samples, GAN could obtain more data. Besides, the transfer learning-related diagnosis model can reuse the previous information in the new task. Consequently, small-sized samples could achieve accurate fault identification.



---

*Source: 1860946-2022-03-30.xml*