# Research on Degradation State Recognition of Planetary Gear Based on Multiscale Information Dimension of SSD and CNN

Authors: Xihui Chen; Liping Peng; Gang Cheng; Chengming Luo
Journal: Complexity
 (2019)

---

## Abstract

Planetary gear is the key part of the transmission system for large complex electromechanical equipment, and in general, a series of degradation states are undergone and evolved into a local fatal fault in its full life cycle. So it is of great significance to recognize the degradation state of planetary gear for the purpose of maintenance repair, predicting development trend, and avoiding sudden fault. This paper proposed a degradation state recognition method of planetary gear based on multiscale information dimension of singular spectrum decomposition (SSD) and convolutional neural network (CNN). SSD can automatically realize the embedding dimension selection and component grouping segmentation, and the original vibration signal being nonlinear and nonstationary can be decomposed into a series of singular spectrum decomposition components (SSDCs), adaptively. Then, the multiscale information dimension which combines multiscale analysis and fractal information dimension is proposed for quantifying and extracting the feature information contained in each SSDC. Finally, CNN is used to achieve the effective recognition of the degradation state of planetary gear. The experimental results show that the proposed method can accurately recognize the degradation state of planetary gear, and the overall recognition rate is up to 97.2%, of which the recognition rate of normal planetary gear reaches 100%.


---

## Body

## 1. Introduction

Planetary gear is the key part of the transmission system for large complex electromechanical equipment, and it usually operates under extreme harsh conditions on the long term. The local fault of planetary gear is easy to produce, and it undergoes a series of degradation states to fatal fault, which directly affects the operational reliability of the electromechanical equipment [1]. Therefore, it is of great significance to accurately recognize the current degradation state of planetary gear for the purpose of maintenance repair, predicting development trend, and avoiding sudden fault.The degradation state recognition of planetary gear is still the problems of pattern classification and fault diagnosis, but it is more difficult than the general fault diagnosis problem [2]. The main reasons are as follows: (1) In actual engineering, the operation condition of planetary gear is relatively harsh. Meanwhile, due to the particularity and complexity of planetary gear structure, its vibration signals have the characteristics of being nonlinear and nonstationary compared with the other transmission mechanisms. (2) The different degradation states of planetary gear still belong to the same fault type, and their degrees are different which results in the fault feature differences being smaller.The effective extraction of the fault feature for planetary gear is the key to realize the degradation state recognition [3]. The traditional fault feature includes time domain features and frequency domain features, but they only have global statistical significance, and they are not suitable to analyze nonlinear and nonstationary signals [4]. With the development of signal processing technology, the method of signal decomposition combining with feature quantification has become the mainstream. The widely used signal decomposition method includes empirical mode decomposition (EMD) [5] and wavelet transform and their improved methods such as ensemble empirical mode decomposition (EEMD) [6], complete ensemble empirical mode decomposition with adaptive noise (CEEMDAN) [7], and dual-tree complex wavelet transform (DTCWT) [8] are also applied to fault feature extraction. But there are still some shortcomings in those methods; the methods developed from EMD can easily fall into modal aliasing when it is used to process the intermittent vibration signal; meanwhile, the obtained decomposition result is as follows: intrinsic mode functions (IMFs) lack the physical meaning. Some methods developed from wavelet transform are restricted to the selection of the wavelet basis functions and decomposition layer. Singular spectrum analysis (SSA) is a kind of nonparametric spectral estimation method based on principal component analysis, and it can capture the high harmonic oscillation shapes based on data adaptive driven, which is suitable for dealing with the vibration signal with nonlinear and nonstationary signals [9, 10]. Traditional SSA is generally divided into four steps: trajectory matrix construction, singular value decomposition (SVD), component grouping, and diagonal averaging. Although some scholars have done some research on SSA, but to determine the appropriate embedding dimension and component grouping criteria is still the key to the effect of SSA. The existing methods need to manually select the window length and embedding dimension, which cannot automatically realize the frequency band division of the decomposed signal, which makes the decomposition result has the shortcomings of physical meaning deficiency and modal aliasing. Therefore, Bonizzi [11] proposed a new adaptive signal processing method: singular spectrum decomposition. This method is a fully adaptive decomposition method, and it can determine the embedding dimension and component grouping criterion based on data driven, and the singular spectrum decomposition components (SSDC) can be reconstructed adaptively from high frequency to low frequency.The original vibration signal with nonlinear and nonstationary characteristics is decomposed into a series of SSDCs by SSD and SSDCs including more feature information that can reflect the degradation state of planetary gear. Therefore, the quantitative extraction of fault features is of vital importance. Analyzing the vibration signal from different time scales can obtain the multidimension expressing of signal information, so multiscale analysis is often applied to signal processing [12]. It can highlight the signal feature in different scales more comprehensively and sufficiently, and it not only reflects the global information but also gives attention to the detail information of vibration signal. The commonly used feature quantization methods include feature frequency, time domain feature, frequency domain feature, entropy [13], and fractal dimension [14]. Fractal dimension can extract the detail signal features, and it reflects the self-similarity of fine structure and statistical significance; it mainly includes Hausdorff dimension, similarity dimension, box dimension, and information dimension [15, 16]. The information dimension can describe the complexity and sparsity of the signal geometric from probabilistic perspective, and it can be used to quantify the fault feature contained in the vibration signal with nonlinear and nonstationary characteristics. Therefore, combining multiscale analysis and information dimension can realize the quantitative extraction of the complexity and sparsity of the vibration signal from different scales.The multiscale information dimension extracted from SSDCs is regarded as fault feature, and the key of the next step is to recognize the degradation state of planetary gear. Traditional pattern recognition methods, such as support vector machine (SVM), back-propagation (BP) neural network, and fuzzy clustering [17, 18], have been applied to the state recognition of mechanical equipment. But those methods still have some shortcomings, which are mainly reflected in the inability to process multidimension data, poor recognition effect with small training samples, and being easy to fall into the local optimum and overfitting. CNN is a multilayer perceptron that applies to recognize two-dimensional feature maps, it is a deep learning network model with multiple hidden layers, and it is beneficial to maintain the relationship among multidimensional data [19]. CNN not only has the advantages of strong fault tolerance, strong self-adaptive ability, and self-learning ability of traditional neural network but also can transfer features layer by layer, and the feature pattern in the training samples can be learned and expressed implicitly. Compared with BP neural network and SVM, CNN can avoid the training falling into the local extreme value, and it has the stronger learning and expression ability to complex features, and it has faster computing speed.The structure of this paper is as follows: In the second section, the mathematical model of degradation state recognition of planetary gear based on multiscale information dimension of SSD and CNN is established. In the third section, the experimental device used in this paper is introduced. In the fourth section, the vibration signals of different degradation states of planetary gear collected by the experimental device are dealt with by the proposed method, which proves the effectiveness of the proposed method for recognizing the degradation state of planetary gear. In the last section, this paper ends with some conclusions.


## 2. Model Building

### 2.1. Singular Spectrum Analysis

SSA is a nonparametric spectrum estimation method to analyze the time domain signal, and it is generally divided into four steps: trajectory matrix construction, SVD, component grouping, and diagonal averaging [20, 21].


#### 2.1.1. Trajectory Matrix Construction

For a nonzero time seriesx1,x2,...,xN with length N, the embedding dimension is m(1<m<N), the one-dimensional time series can be converted into multidimensional signal space X~1T,X~2T,...,X~KT, and the obtained lag vector can be expressed as X~i=(xi,xi+1,...,xi+m-1). The obtained trajectory matrix can be expressed as follows:(1)X=X~1T,X~2T,…,X~KT=x1x2⋯xKx2x3⋯xK+1⋮⋮⋯⋮xmxm+1⋯xN,where the obtained trajectory matrix X is a Hankel matrix, and it contains the information with trend, oscillation, and noise. It can be seen that the choice of embedding dimension is an important parameter for the trajectory matrix construction, which affects the effect of SVD.


#### 2.1.2. Singular Value Decomposition (SVD)

For the obtained trajectory matrixX(m×k), SVD is carried out. The SVD of trajectory matrix X can be expressed as follows: X=UDVT, where U and V are normalization matrix, and D is the diagonal matrix constructed by the obtained singular value σi. According to the theory of SVD, the singular value σi can be obtained by eigenvalue decomposition of the covariance matrix of trajectory matrix X. The covariance matrix of trajectory matrix X can be expressed as C=(1/K)XXT, where K is the dimension. The eigenvalues of the covariance matrix C can be obtained and sorted from large to small, and they are λ=λ1,λ2,...,λK; further, the obtained singular values can be expressed as σi=Kλi. The matrix X can be decomposed into a sum of submatrix with order 1, and X=∑i=1MXi=∑i=1Mσiuivi, where Xi represents thei-th matrix obtained by SVD.


#### 2.1.3. Component Grouping

The component grouping process is to divide matrix setX1,X2,...,XM into S group disjoint matrix subset XI1,XI2,...,XIS, and the sum of each group of matrix XIs can be calculated, if Is=i1,i2,...,ip; then, XIs=Xi1+Xi2+...+Xip. The above process is also called as eigentriple grouping, it is the key step of SSA, and the components containing different signals (trend, oscillation, noise, etc.) are decomposed into different groups. But the determination of component grouping criterion is generally based on artificial setting.


#### 2.1.4. Diagonal Averaging

This process is the last step of SSA, diagonal averaging is carried out for the new grouping matrixXIS, and the final reconstruction components that are the same as the data points of the original signal can be obtained. This process is called diagonal averaging or Hankelization. The corresponding relationship between the final reconstruction components Y(i) by diagonal averaging and the each data points of new grouping matrix XIS is as follows:(2)y~n=1n∑j=1naj,n-j+1,1≤n≤m,1L∑j=1Laj,n-j+1,m≤n≤K,1N-n+1∑j=n-K+1Laj,n-j+1,K≤n≤N,where y~n represents the i-th SSDC, aj,n represents the data point of new grouping matrix XIs obtained by SVD and component grouping, and K=N-m+1.The above equation is the main process of SSA; it can be seen that the critical parameters that determine the decomposition components quality of SSA are the selection of embedding dimensionm and component grouping criterion.


### 2.2. Singular Spectrum Decomposition

In order to overcome the difficulty of SSA in selecting the embedding dimension and how to carry out the component grouping, a new adaptive signal processing method named SSD is proposed. It can automatically select the embedding dimension in the iterative process, and the frequency band of the obtained SSDCs can be segmented actively by automatic component grouping. It is a completely data driven decomposition method [11, 22].


#### 2.2.1. Improved Construction of Trajectory Matrix

Unlike (1), a new improved method is used to construct the trajectory matrix. The improved construction process of trajectory matrix is as follows: For a nonzero time series x1,x2,...,xN with length N, the embedding dimension is m(1<m<N), and the one-dimensional time series can be converted into multidimensional signal space X=X^1;X^2;...;X^m, where X^i=xi,xi+1,...,xN,x1,x2,...,xi-1,i=1,2,...,m. Assuming time series is 1,2,3,4,5,6, and the embedding dimension m is set as 4, then, the obtained trajectory matrix is as follows:(3)X=123456234561345612456123.It can be seen that the left four columns of (3) corresponds to the trajectory matrix of the standard SSA. The trajectory matrix of SSD is improved by comparing SSA, the main purpose is to enhance the oscillation component in original signal, and the energy of residual signal can be reduced to display the properties of useful signal.


#### 2.2.2. Adaptive Selection of Embedding Dimension

The selection of embedding dimension is important for constructing the trajectory matrix of SSD, and the adaptive selection criterion of embedding dimension is formulated. It can select the embedding dimension adaptively in each iteration process to calculate the final SSDDs based on data driven. The adaptive selection criterion of embedding dimension is expressed as follows:(1) According to SSD theory [22], a SSDC can be obtained through one iteration calculation. Assuming the residual component in j-th iteration is vj(n)=x(n)-∑k=1j-1vk(n), v0(n)=x(n), and its power spectrum density (PSD) is calculated. The frequency fmax corresponding to maximum peak in PSD can be obtained.(2) In the first iteration, namely,j=1, if the value of fmax is smaller (this is measured by whether fmax/fs is less than 0.01, where fs is the sampling frequency), it is shown that the residual signal is considered as a trend signal, and the embedding dimension m is set as N/3.(3) In other cases, if it is not the first iteration, namely,j>1, then m is set as l∙(Fs/fmax), where l is the ratio factor for adjusting the average period of the desired signal and the window length; in general, l is set as 1.2.


#### 2.2.3. Grouping and Reconstruction of thej-th SSDC

SSDCs are reconstructed from high frequency to low frequency by a series of iterations in SSD process. In the first iteration (j=1), if the detected component is a trend signal, the first left and right eigenvectors are only used to reconstruct the SSDC g(1)(n), and g(1)(n)=dia(X1)=dia(σ1u1v1T), where dia represents the calculation process of diagonal averaging. In addition, for the j-th iteration, the SSDC g(j)(n) must be able to describe a time scale with clear physical meaning. In this sense, the following rules are defined: The main frequency components of the analyzed signal are concentrated in the frequency band fmax-Δf,fmax+Δf, where Δf represents the half band width of the main peak in the PSD of the residual signal. Therefore, a subset Ij=(Ij=i1,i2,...,ip) is created from all eigenvalues set, and the determining principle of subset Ij is as follows: All eigenvalues which correspond to the left eigenvectors have the prominent dominant frequency in frequency band fmax-Δf,fmax+Δf and eigenvalues with the greatest contribution to the main peak energy of the analyzed signal are selected as a group; then, the SSDC corresponding to subset Ij can be reconstructed according to the diagonal averaging of the matrix XIj=Xi1+Xi2+⋯+Xip. The half band width of the main peak Δf in PSD is related to the average time span of the oscillating signal of the SSDC. In order to better estimate Δf, a spectrum model with the superposed Gaussian function is constructed to describe the distribution of PSD. This model is defined as the sum of three Gaussian functions, and each function represents a spectrum peak:(4)rf,θ=∑i=13Aie-f-ui2/2σi2,where Ai is the amplitude of the i -th Gaussian function, ui is its location, and σi is its width; θ=AσT is a parameter vector, and it satisfies A=[A1,A2,A3] and σ=[σ1,σ2,σ3]. The first Gaussian function closes to the frequency fmax corresponding to the main spectrum peak, and the second Gaussian function closes to the frequency f2 corresponding to the subspectrum peak, and the third Gaussian function closes to the frequency corresponding to any peak between main spectrum peak and subspectrum peak. Therefore, the following can be obtained:(5)u1=fmax,u2=f2,u3=f3=fmax+f22.The model parameterAi can be obtained by the weighted least square method, and the initial parameter values of the model can be set as follows:(6)A10=12PSDfmax,σ10=f:PSDf=23PSDfmaxA20=12PSDf2,σ20=f:PSDf=23PSDf2A30=14PSDf3,σ30=4fmax-f2.The optimal value of model parameterAi is determined by Levenberg-Marquardt method. The estimated value of σ1 is given, and the bandwidth Δf=2.5σ1 of the main spectrum peak can be obtained. Choosing XIj by the above method can determine the main eigenvalues with respect to the noise effects. In the j-th iteration, the signal components with different scales mismatching to the frequency band fmax-Δf,fmax+Δf in this iteration are automatically discarded, and they can be decomposed in subsequent iteration. In further, in order to reconstruct the j-th SSDC, thej-th iteration begins. The scale factor a^ is used to adjust the difference between g(1)(n) and residual signal v(j)(n), and it is as follows:(7)a^=mina⁡vjn-ag~1n22,wherea^=gTvgTgandg~jn=a^gjn.


#### 2.2.4. Stopping Criterion for Iteration

The SSDCsg~(j)(n) obtained by iterations are separated from original vibration signal, and the residual signal is v(j+1)(n)=v(j)(n)-g~(j)(n). The normalized mean square error (NMSE) is calculated between the residual signal and original vibration signal, and it is as follows:(8)NMSEj=∑i=1Nvj+1i2∑i=1Nxi2.NMSE is less than a certain threshold value set as the stopping criterion for iteration, and this threshold value can be set as 1%. If NSME is less than the threshold value, then the iteration decomposition is terminated. If NSME is greater than the threshold value, the residual signal is regarded as the original signal to repeat the above iteration decomposition process until the stopping criterion for iteration is satisfied. After satisfying the stopping criterion for iteration, the final decomposed SSDCs are as follows:(9)xn=∑k=1Mg~kn+vM+1n,where M is the number of the decomposed SSDCs and v(M+1)(n) is the residual component when iteration decomposition is terminated.


### 2.3. Multiscale Fractal Information Dimension

The original vibration signal with nonlinear and nonstationary characteristics is decomposed into a series of SSDCs by SSD, and combining multiscale analysis and fractal information dimension can realize the quantitative extraction of the complexity and sparsity of the vibration signal from different scales [12].


#### 2.3.1. Multiscale Analysis

Multiscale analysis is that a new data point can be obtained by average processing of the adjacentτ data points, and the new data points can form a new time domain signal with τ scale. For a time domain signal x1,x2,...,xN, the value of scale factor τ is set, and a new time domain signal yjτ can be obtained in τ scale, and the specific process is shown as (10) and Figure 1. By changing the scale factor τ, the expression of the original signal in different scales can be obtained:(10)yjτ=1τ∑h=j-1τ+1jτxh1≤j≤Hτ.Figure 1
Multiscale analysis.


#### 2.3.2. Fractal Information Dimension

Fractal information dimension describes the complexity and sparsity of the signal geometric from probabilistic perspective. The specific calculation process of information dimension is as follows: For the signal sequencey=y1,y2,…,yN, in order to reduce the influence of noise on the signal, the difference of adjacent data points in the effective signal length is used as the reconstruction signal, that is, s0=si+1-si,i=1,2,…,N-1. Assuming that the space of the reconstruction signal is filled by a series of boxes with a length of δ, and all boxes are numbered in the calculation process. If the probability that the data points of the reconstruction signal fall into the i-th box is Pi, then, the information entropy can be expressed as(11)S=∑i=1N-1s0ipi=s0iSHI=-∑i=1N-1pilg⁡pi.If the information dimension satisfiesHI~lg⁡δDI, the information dimension of X can be defined as follows:(12)DI=-limδ→0⁡HIδlg⁡δ.


### 2.4. Convolutional Neural Network

CNN is generally composed of input layer, convolution layer, pooling layer, full connection layer, and classifier layer, and its basic structure is shown as Figure2 [23–25].Figure 2
The basic structure of CNN.


#### 2.4.1. Convolution Layer

The convolution layer contains a set of convolution kernels with equal size learned by data driven. For different feature density, the convolution kernels with the fixed size convolution step are convoluted with the input feature matrix, and the convolution feature map can be formed by nonlinear activation function and bias, which represents the response to input feature, and the calculation process of the convolution feature map can be expressed as follows:(13)anl=f∑∀maml-1∗km,nl+bnl,where anl and aml-1 are the n-th feature map and m-th feature map of the l-th layer and l+1-th layer, respectively. km,nl is a convolution kernel between two feature maps, bnl is the bias, and f is the nonlinear activation function, and in here, Sigmoid function is used and shown as follows:(14)fx=11+e-x.


#### 2.4.2. Pooling Layer

The pooling layer is usually cascaded with convolution layer, and its function is to reduce the dimension of convolution feature map. The pooling feature map is formed by the downsampling of convolution feature map, and the calculation process of pooling feature map can be expressed as follows:(15)anl=fknl×1s2∑s×sanl-1+bnl,where s is the template size of downsampling and kml is the template weight. According to different sampling method, the pooling computing methods include maximum pooling, average pooling, and random pooling. Pooling is the aggregation statistics of the feature in a continuous region, and the regional feature is represented by the maximum value and average value and so on. In this paper, the maximum pooling is used.


#### 2.4.3. Full Connection Layer

The full connection layer adopts the full connection mode, and the vector transformation of pooling feature map is processed. The two-dimensional feature matrix of this layer is stretched into one-dimension feature vector, which is convenient for the subsequent output layer calculation. The full connection layer is still equivalent to the convolution layer, and the difference is that the arrangement of the output is different, or the elements of pooling feature map are arranged directly into one-dimension vector.


#### 2.4.4. Classifier Layer

The classifier layer mainly implements two functions: (1) In training stage of CNN, it provides the labels of the training data, and the training error is calculated. The weights are adjusted by gradient descent and reverse propagation, and the training error is reduced to minimum which completes the training process of CNN. (2) In testing stage of CNN, the output probability of testing samples is calculated, and the state of testing samples is determined according to the criterion of maximum probability. At present, the most commonly used classification method in classifier layer is logistic method and Softmax method. Because the degradation state recognition of planetary gear is a multipattern classification problem, in this paper, the Softmax method is used.


## 2.1. Singular Spectrum Analysis

SSA is a nonparametric spectrum estimation method to analyze the time domain signal, and it is generally divided into four steps: trajectory matrix construction, SVD, component grouping, and diagonal averaging [20, 21].


### 2.1.1. Trajectory Matrix Construction

For a nonzero time seriesx1,x2,...,xN with length N, the embedding dimension is m(1<m<N), the one-dimensional time series can be converted into multidimensional signal space X~1T,X~2T,...,X~KT, and the obtained lag vector can be expressed as X~i=(xi,xi+1,...,xi+m-1). The obtained trajectory matrix can be expressed as follows:(1)X=X~1T,X~2T,…,X~KT=x1x2⋯xKx2x3⋯xK+1⋮⋮⋯⋮xmxm+1⋯xN,where the obtained trajectory matrix X is a Hankel matrix, and it contains the information with trend, oscillation, and noise. It can be seen that the choice of embedding dimension is an important parameter for the trajectory matrix construction, which affects the effect of SVD.


### 2.1.2. Singular Value Decomposition (SVD)

For the obtained trajectory matrixX(m×k), SVD is carried out. The SVD of trajectory matrix X can be expressed as follows: X=UDVT, where U and V are normalization matrix, and D is the diagonal matrix constructed by the obtained singular value σi. According to the theory of SVD, the singular value σi can be obtained by eigenvalue decomposition of the covariance matrix of trajectory matrix X. The covariance matrix of trajectory matrix X can be expressed as C=(1/K)XXT, where K is the dimension. The eigenvalues of the covariance matrix C can be obtained and sorted from large to small, and they are λ=λ1,λ2,...,λK; further, the obtained singular values can be expressed as σi=Kλi. The matrix X can be decomposed into a sum of submatrix with order 1, and X=∑i=1MXi=∑i=1Mσiuivi, where Xi represents thei-th matrix obtained by SVD.


### 2.1.3. Component Grouping

The component grouping process is to divide matrix setX1,X2,...,XM into S group disjoint matrix subset XI1,XI2,...,XIS, and the sum of each group of matrix XIs can be calculated, if Is=i1,i2,...,ip; then, XIs=Xi1+Xi2+...+Xip. The above process is also called as eigentriple grouping, it is the key step of SSA, and the components containing different signals (trend, oscillation, noise, etc.) are decomposed into different groups. But the determination of component grouping criterion is generally based on artificial setting.


### 2.1.4. Diagonal Averaging

This process is the last step of SSA, diagonal averaging is carried out for the new grouping matrixXIS, and the final reconstruction components that are the same as the data points of the original signal can be obtained. This process is called diagonal averaging or Hankelization. The corresponding relationship between the final reconstruction components Y(i) by diagonal averaging and the each data points of new grouping matrix XIS is as follows:(2)y~n=1n∑j=1naj,n-j+1,1≤n≤m,1L∑j=1Laj,n-j+1,m≤n≤K,1N-n+1∑j=n-K+1Laj,n-j+1,K≤n≤N,where y~n represents the i-th SSDC, aj,n represents the data point of new grouping matrix XIs obtained by SVD and component grouping, and K=N-m+1.The above equation is the main process of SSA; it can be seen that the critical parameters that determine the decomposition components quality of SSA are the selection of embedding dimensionm and component grouping criterion.


## 2.1.1. Trajectory Matrix Construction

For a nonzero time seriesx1,x2,...,xN with length N, the embedding dimension is m(1<m<N), the one-dimensional time series can be converted into multidimensional signal space X~1T,X~2T,...,X~KT, and the obtained lag vector can be expressed as X~i=(xi,xi+1,...,xi+m-1). The obtained trajectory matrix can be expressed as follows:(1)X=X~1T,X~2T,…,X~KT=x1x2⋯xKx2x3⋯xK+1⋮⋮⋯⋮xmxm+1⋯xN,where the obtained trajectory matrix X is a Hankel matrix, and it contains the information with trend, oscillation, and noise. It can be seen that the choice of embedding dimension is an important parameter for the trajectory matrix construction, which affects the effect of SVD.


## 2.1.2. Singular Value Decomposition (SVD)

For the obtained trajectory matrixX(m×k), SVD is carried out. The SVD of trajectory matrix X can be expressed as follows: X=UDVT, where U and V are normalization matrix, and D is the diagonal matrix constructed by the obtained singular value σi. According to the theory of SVD, the singular value σi can be obtained by eigenvalue decomposition of the covariance matrix of trajectory matrix X. The covariance matrix of trajectory matrix X can be expressed as C=(1/K)XXT, where K is the dimension. The eigenvalues of the covariance matrix C can be obtained and sorted from large to small, and they are λ=λ1,λ2,...,λK; further, the obtained singular values can be expressed as σi=Kλi. The matrix X can be decomposed into a sum of submatrix with order 1, and X=∑i=1MXi=∑i=1Mσiuivi, where Xi represents thei-th matrix obtained by SVD.


## 2.1.3. Component Grouping

The component grouping process is to divide matrix setX1,X2,...,XM into S group disjoint matrix subset XI1,XI2,...,XIS, and the sum of each group of matrix XIs can be calculated, if Is=i1,i2,...,ip; then, XIs=Xi1+Xi2+...+Xip. The above process is also called as eigentriple grouping, it is the key step of SSA, and the components containing different signals (trend, oscillation, noise, etc.) are decomposed into different groups. But the determination of component grouping criterion is generally based on artificial setting.


## 2.1.4. Diagonal Averaging

This process is the last step of SSA, diagonal averaging is carried out for the new grouping matrixXIS, and the final reconstruction components that are the same as the data points of the original signal can be obtained. This process is called diagonal averaging or Hankelization. The corresponding relationship between the final reconstruction components Y(i) by diagonal averaging and the each data points of new grouping matrix XIS is as follows:(2)y~n=1n∑j=1naj,n-j+1,1≤n≤m,1L∑j=1Laj,n-j+1,m≤n≤K,1N-n+1∑j=n-K+1Laj,n-j+1,K≤n≤N,where y~n represents the i-th SSDC, aj,n represents the data point of new grouping matrix XIs obtained by SVD and component grouping, and K=N-m+1.The above equation is the main process of SSA; it can be seen that the critical parameters that determine the decomposition components quality of SSA are the selection of embedding dimensionm and component grouping criterion.


## 2.2. Singular Spectrum Decomposition

In order to overcome the difficulty of SSA in selecting the embedding dimension and how to carry out the component grouping, a new adaptive signal processing method named SSD is proposed. It can automatically select the embedding dimension in the iterative process, and the frequency band of the obtained SSDCs can be segmented actively by automatic component grouping. It is a completely data driven decomposition method [11, 22].


### 2.2.1. Improved Construction of Trajectory Matrix

Unlike (1), a new improved method is used to construct the trajectory matrix. The improved construction process of trajectory matrix is as follows: For a nonzero time series x1,x2,...,xN with length N, the embedding dimension is m(1<m<N), and the one-dimensional time series can be converted into multidimensional signal space X=X^1;X^2;...;X^m, where X^i=xi,xi+1,...,xN,x1,x2,...,xi-1,i=1,2,...,m. Assuming time series is 1,2,3,4,5,6, and the embedding dimension m is set as 4, then, the obtained trajectory matrix is as follows:(3)X=123456234561345612456123.It can be seen that the left four columns of (3) corresponds to the trajectory matrix of the standard SSA. The trajectory matrix of SSD is improved by comparing SSA, the main purpose is to enhance the oscillation component in original signal, and the energy of residual signal can be reduced to display the properties of useful signal.


### 2.2.2. Adaptive Selection of Embedding Dimension

The selection of embedding dimension is important for constructing the trajectory matrix of SSD, and the adaptive selection criterion of embedding dimension is formulated. It can select the embedding dimension adaptively in each iteration process to calculate the final SSDDs based on data driven. The adaptive selection criterion of embedding dimension is expressed as follows:(1) According to SSD theory [22], a SSDC can be obtained through one iteration calculation. Assuming the residual component in j-th iteration is vj(n)=x(n)-∑k=1j-1vk(n), v0(n)=x(n), and its power spectrum density (PSD) is calculated. The frequency fmax corresponding to maximum peak in PSD can be obtained.(2) In the first iteration, namely,j=1, if the value of fmax is smaller (this is measured by whether fmax/fs is less than 0.01, where fs is the sampling frequency), it is shown that the residual signal is considered as a trend signal, and the embedding dimension m is set as N/3.(3) In other cases, if it is not the first iteration, namely,j>1, then m is set as l∙(Fs/fmax), where l is the ratio factor for adjusting the average period of the desired signal and the window length; in general, l is set as 1.2.


### 2.2.3. Grouping and Reconstruction of thej-th SSDC

SSDCs are reconstructed from high frequency to low frequency by a series of iterations in SSD process. In the first iteration (j=1), if the detected component is a trend signal, the first left and right eigenvectors are only used to reconstruct the SSDC g(1)(n), and g(1)(n)=dia(X1)=dia(σ1u1v1T), where dia represents the calculation process of diagonal averaging. In addition, for the j-th iteration, the SSDC g(j)(n) must be able to describe a time scale with clear physical meaning. In this sense, the following rules are defined: The main frequency components of the analyzed signal are concentrated in the frequency band fmax-Δf,fmax+Δf, where Δf represents the half band width of the main peak in the PSD of the residual signal. Therefore, a subset Ij=(Ij=i1,i2,...,ip) is created from all eigenvalues set, and the determining principle of subset Ij is as follows: All eigenvalues which correspond to the left eigenvectors have the prominent dominant frequency in frequency band fmax-Δf,fmax+Δf and eigenvalues with the greatest contribution to the main peak energy of the analyzed signal are selected as a group; then, the SSDC corresponding to subset Ij can be reconstructed according to the diagonal averaging of the matrix XIj=Xi1+Xi2+⋯+Xip. The half band width of the main peak Δf in PSD is related to the average time span of the oscillating signal of the SSDC. In order to better estimate Δf, a spectrum model with the superposed Gaussian function is constructed to describe the distribution of PSD. This model is defined as the sum of three Gaussian functions, and each function represents a spectrum peak:(4)rf,θ=∑i=13Aie-f-ui2/2σi2,where Ai is the amplitude of the i -th Gaussian function, ui is its location, and σi is its width; θ=AσT is a parameter vector, and it satisfies A=[A1,A2,A3] and σ=[σ1,σ2,σ3]. The first Gaussian function closes to the frequency fmax corresponding to the main spectrum peak, and the second Gaussian function closes to the frequency f2 corresponding to the subspectrum peak, and the third Gaussian function closes to the frequency corresponding to any peak between main spectrum peak and subspectrum peak. Therefore, the following can be obtained:(5)u1=fmax,u2=f2,u3=f3=fmax+f22.The model parameterAi can be obtained by the weighted least square method, and the initial parameter values of the model can be set as follows:(6)A10=12PSDfmax,σ10=f:PSDf=23PSDfmaxA20=12PSDf2,σ20=f:PSDf=23PSDf2A30=14PSDf3,σ30=4fmax-f2.The optimal value of model parameterAi is determined by Levenberg-Marquardt method. The estimated value of σ1 is given, and the bandwidth Δf=2.5σ1 of the main spectrum peak can be obtained. Choosing XIj by the above method can determine the main eigenvalues with respect to the noise effects. In the j-th iteration, the signal components with different scales mismatching to the frequency band fmax-Δf,fmax+Δf in this iteration are automatically discarded, and they can be decomposed in subsequent iteration. In further, in order to reconstruct the j-th SSDC, thej-th iteration begins. The scale factor a^ is used to adjust the difference between g(1)(n) and residual signal v(j)(n), and it is as follows:(7)a^=mina⁡vjn-ag~1n22,wherea^=gTvgTgandg~jn=a^gjn.


### 2.2.4. Stopping Criterion for Iteration

The SSDCsg~(j)(n) obtained by iterations are separated from original vibration signal, and the residual signal is v(j+1)(n)=v(j)(n)-g~(j)(n). The normalized mean square error (NMSE) is calculated between the residual signal and original vibration signal, and it is as follows:(8)NMSEj=∑i=1Nvj+1i2∑i=1Nxi2.NMSE is less than a certain threshold value set as the stopping criterion for iteration, and this threshold value can be set as 1%. If NSME is less than the threshold value, then the iteration decomposition is terminated. If NSME is greater than the threshold value, the residual signal is regarded as the original signal to repeat the above iteration decomposition process until the stopping criterion for iteration is satisfied. After satisfying the stopping criterion for iteration, the final decomposed SSDCs are as follows:(9)xn=∑k=1Mg~kn+vM+1n,where M is the number of the decomposed SSDCs and v(M+1)(n) is the residual component when iteration decomposition is terminated.


## 2.2.1. Improved Construction of Trajectory Matrix

Unlike (1), a new improved method is used to construct the trajectory matrix. The improved construction process of trajectory matrix is as follows: For a nonzero time series x1,x2,...,xN with length N, the embedding dimension is m(1<m<N), and the one-dimensional time series can be converted into multidimensional signal space X=X^1;X^2;...;X^m, where X^i=xi,xi+1,...,xN,x1,x2,...,xi-1,i=1,2,...,m. Assuming time series is 1,2,3,4,5,6, and the embedding dimension m is set as 4, then, the obtained trajectory matrix is as follows:(3)X=123456234561345612456123.It can be seen that the left four columns of (3) corresponds to the trajectory matrix of the standard SSA. The trajectory matrix of SSD is improved by comparing SSA, the main purpose is to enhance the oscillation component in original signal, and the energy of residual signal can be reduced to display the properties of useful signal.


## 2.2.2. Adaptive Selection of Embedding Dimension

The selection of embedding dimension is important for constructing the trajectory matrix of SSD, and the adaptive selection criterion of embedding dimension is formulated. It can select the embedding dimension adaptively in each iteration process to calculate the final SSDDs based on data driven. The adaptive selection criterion of embedding dimension is expressed as follows:(1) According to SSD theory [22], a SSDC can be obtained through one iteration calculation. Assuming the residual component in j-th iteration is vj(n)=x(n)-∑k=1j-1vk(n), v0(n)=x(n), and its power spectrum density (PSD) is calculated. The frequency fmax corresponding to maximum peak in PSD can be obtained.(2) In the first iteration, namely,j=1, if the value of fmax is smaller (this is measured by whether fmax/fs is less than 0.01, where fs is the sampling frequency), it is shown that the residual signal is considered as a trend signal, and the embedding dimension m is set as N/3.(3) In other cases, if it is not the first iteration, namely,j>1, then m is set as l∙(Fs/fmax), where l is the ratio factor for adjusting the average period of the desired signal and the window length; in general, l is set as 1.2.


## 2.2.3. Grouping and Reconstruction of thej-th SSDC

SSDCs are reconstructed from high frequency to low frequency by a series of iterations in SSD process. In the first iteration (j=1), if the detected component is a trend signal, the first left and right eigenvectors are only used to reconstruct the SSDC g(1)(n), and g(1)(n)=dia(X1)=dia(σ1u1v1T), where dia represents the calculation process of diagonal averaging. In addition, for the j-th iteration, the SSDC g(j)(n) must be able to describe a time scale with clear physical meaning. In this sense, the following rules are defined: The main frequency components of the analyzed signal are concentrated in the frequency band fmax-Δf,fmax+Δf, where Δf represents the half band width of the main peak in the PSD of the residual signal. Therefore, a subset Ij=(Ij=i1,i2,...,ip) is created from all eigenvalues set, and the determining principle of subset Ij is as follows: All eigenvalues which correspond to the left eigenvectors have the prominent dominant frequency in frequency band fmax-Δf,fmax+Δf and eigenvalues with the greatest contribution to the main peak energy of the analyzed signal are selected as a group; then, the SSDC corresponding to subset Ij can be reconstructed according to the diagonal averaging of the matrix XIj=Xi1+Xi2+⋯+Xip. The half band width of the main peak Δf in PSD is related to the average time span of the oscillating signal of the SSDC. In order to better estimate Δf, a spectrum model with the superposed Gaussian function is constructed to describe the distribution of PSD. This model is defined as the sum of three Gaussian functions, and each function represents a spectrum peak:(4)rf,θ=∑i=13Aie-f-ui2/2σi2,where Ai is the amplitude of the i -th Gaussian function, ui is its location, and σi is its width; θ=AσT is a parameter vector, and it satisfies A=[A1,A2,A3] and σ=[σ1,σ2,σ3]. The first Gaussian function closes to the frequency fmax corresponding to the main spectrum peak, and the second Gaussian function closes to the frequency f2 corresponding to the subspectrum peak, and the third Gaussian function closes to the frequency corresponding to any peak between main spectrum peak and subspectrum peak. Therefore, the following can be obtained:(5)u1=fmax,u2=f2,u3=f3=fmax+f22.The model parameterAi can be obtained by the weighted least square method, and the initial parameter values of the model can be set as follows:(6)A10=12PSDfmax,σ10=f:PSDf=23PSDfmaxA20=12PSDf2,σ20=f:PSDf=23PSDf2A30=14PSDf3,σ30=4fmax-f2.The optimal value of model parameterAi is determined by Levenberg-Marquardt method. The estimated value of σ1 is given, and the bandwidth Δf=2.5σ1 of the main spectrum peak can be obtained. Choosing XIj by the above method can determine the main eigenvalues with respect to the noise effects. In the j-th iteration, the signal components with different scales mismatching to the frequency band fmax-Δf,fmax+Δf in this iteration are automatically discarded, and they can be decomposed in subsequent iteration. In further, in order to reconstruct the j-th SSDC, thej-th iteration begins. The scale factor a^ is used to adjust the difference between g(1)(n) and residual signal v(j)(n), and it is as follows:(7)a^=mina⁡vjn-ag~1n22,wherea^=gTvgTgandg~jn=a^gjn.


## 2.2.4. Stopping Criterion for Iteration

The SSDCsg~(j)(n) obtained by iterations are separated from original vibration signal, and the residual signal is v(j+1)(n)=v(j)(n)-g~(j)(n). The normalized mean square error (NMSE) is calculated between the residual signal and original vibration signal, and it is as follows:(8)NMSEj=∑i=1Nvj+1i2∑i=1Nxi2.NMSE is less than a certain threshold value set as the stopping criterion for iteration, and this threshold value can be set as 1%. If NSME is less than the threshold value, then the iteration decomposition is terminated. If NSME is greater than the threshold value, the residual signal is regarded as the original signal to repeat the above iteration decomposition process until the stopping criterion for iteration is satisfied. After satisfying the stopping criterion for iteration, the final decomposed SSDCs are as follows:(9)xn=∑k=1Mg~kn+vM+1n,where M is the number of the decomposed SSDCs and v(M+1)(n) is the residual component when iteration decomposition is terminated.


## 2.3. Multiscale Fractal Information Dimension

The original vibration signal with nonlinear and nonstationary characteristics is decomposed into a series of SSDCs by SSD, and combining multiscale analysis and fractal information dimension can realize the quantitative extraction of the complexity and sparsity of the vibration signal from different scales [12].


### 2.3.1. Multiscale Analysis

Multiscale analysis is that a new data point can be obtained by average processing of the adjacentτ data points, and the new data points can form a new time domain signal with τ scale. For a time domain signal x1,x2,...,xN, the value of scale factor τ is set, and a new time domain signal yjτ can be obtained in τ scale, and the specific process is shown as (10) and Figure 1. By changing the scale factor τ, the expression of the original signal in different scales can be obtained:(10)yjτ=1τ∑h=j-1τ+1jτxh1≤j≤Hτ.Figure 1
Multiscale analysis.


### 2.3.2. Fractal Information Dimension

Fractal information dimension describes the complexity and sparsity of the signal geometric from probabilistic perspective. The specific calculation process of information dimension is as follows: For the signal sequencey=y1,y2,…,yN, in order to reduce the influence of noise on the signal, the difference of adjacent data points in the effective signal length is used as the reconstruction signal, that is, s0=si+1-si,i=1,2,…,N-1. Assuming that the space of the reconstruction signal is filled by a series of boxes with a length of δ, and all boxes are numbered in the calculation process. If the probability that the data points of the reconstruction signal fall into the i-th box is Pi, then, the information entropy can be expressed as(11)S=∑i=1N-1s0ipi=s0iSHI=-∑i=1N-1pilg⁡pi.If the information dimension satisfiesHI~lg⁡δDI, the information dimension of X can be defined as follows:(12)DI=-limδ→0⁡HIδlg⁡δ.


## 2.3.1. Multiscale Analysis

Multiscale analysis is that a new data point can be obtained by average processing of the adjacentτ data points, and the new data points can form a new time domain signal with τ scale. For a time domain signal x1,x2,...,xN, the value of scale factor τ is set, and a new time domain signal yjτ can be obtained in τ scale, and the specific process is shown as (10) and Figure 1. By changing the scale factor τ, the expression of the original signal in different scales can be obtained:(10)yjτ=1τ∑h=j-1τ+1jτxh1≤j≤Hτ.Figure 1
Multiscale analysis.


## 2.3.2. Fractal Information Dimension

Fractal information dimension describes the complexity and sparsity of the signal geometric from probabilistic perspective. The specific calculation process of information dimension is as follows: For the signal sequencey=y1,y2,…,yN, in order to reduce the influence of noise on the signal, the difference of adjacent data points in the effective signal length is used as the reconstruction signal, that is, s0=si+1-si,i=1,2,…,N-1. Assuming that the space of the reconstruction signal is filled by a series of boxes with a length of δ, and all boxes are numbered in the calculation process. If the probability that the data points of the reconstruction signal fall into the i-th box is Pi, then, the information entropy can be expressed as(11)S=∑i=1N-1s0ipi=s0iSHI=-∑i=1N-1pilg⁡pi.If the information dimension satisfiesHI~lg⁡δDI, the information dimension of X can be defined as follows:(12)DI=-limδ→0⁡HIδlg⁡δ.


## 2.4. Convolutional Neural Network

CNN is generally composed of input layer, convolution layer, pooling layer, full connection layer, and classifier layer, and its basic structure is shown as Figure2 [23–25].Figure 2
The basic structure of CNN.


### 2.4.1. Convolution Layer

The convolution layer contains a set of convolution kernels with equal size learned by data driven. For different feature density, the convolution kernels with the fixed size convolution step are convoluted with the input feature matrix, and the convolution feature map can be formed by nonlinear activation function and bias, which represents the response to input feature, and the calculation process of the convolution feature map can be expressed as follows:(13)anl=f∑∀maml-1∗km,nl+bnl,where anl and aml-1 are the n-th feature map and m-th feature map of the l-th layer and l+1-th layer, respectively. km,nl is a convolution kernel between two feature maps, bnl is the bias, and f is the nonlinear activation function, and in here, Sigmoid function is used and shown as follows:(14)fx=11+e-x.


### 2.4.2. Pooling Layer

The pooling layer is usually cascaded with convolution layer, and its function is to reduce the dimension of convolution feature map. The pooling feature map is formed by the downsampling of convolution feature map, and the calculation process of pooling feature map can be expressed as follows:(15)anl=fknl×1s2∑s×sanl-1+bnl,where s is the template size of downsampling and kml is the template weight. According to different sampling method, the pooling computing methods include maximum pooling, average pooling, and random pooling. Pooling is the aggregation statistics of the feature in a continuous region, and the regional feature is represented by the maximum value and average value and so on. In this paper, the maximum pooling is used.


### 2.4.3. Full Connection Layer

The full connection layer adopts the full connection mode, and the vector transformation of pooling feature map is processed. The two-dimensional feature matrix of this layer is stretched into one-dimension feature vector, which is convenient for the subsequent output layer calculation. The full connection layer is still equivalent to the convolution layer, and the difference is that the arrangement of the output is different, or the elements of pooling feature map are arranged directly into one-dimension vector.


### 2.4.4. Classifier Layer

The classifier layer mainly implements two functions: (1) In training stage of CNN, it provides the labels of the training data, and the training error is calculated. The weights are adjusted by gradient descent and reverse propagation, and the training error is reduced to minimum which completes the training process of CNN. (2) In testing stage of CNN, the output probability of testing samples is calculated, and the state of testing samples is determined according to the criterion of maximum probability. At present, the most commonly used classification method in classifier layer is logistic method and Softmax method. Because the degradation state recognition of planetary gear is a multipattern classification problem, in this paper, the Softmax method is used.


## 2.4.1. Convolution Layer

The convolution layer contains a set of convolution kernels with equal size learned by data driven. For different feature density, the convolution kernels with the fixed size convolution step are convoluted with the input feature matrix, and the convolution feature map can be formed by nonlinear activation function and bias, which represents the response to input feature, and the calculation process of the convolution feature map can be expressed as follows:(13)anl=f∑∀maml-1∗km,nl+bnl,where anl and aml-1 are the n-th feature map and m-th feature map of the l-th layer and l+1-th layer, respectively. km,nl is a convolution kernel between two feature maps, bnl is the bias, and f is the nonlinear activation function, and in here, Sigmoid function is used and shown as follows:(14)fx=11+e-x.


## 2.4.2. Pooling Layer

The pooling layer is usually cascaded with convolution layer, and its function is to reduce the dimension of convolution feature map. The pooling feature map is formed by the downsampling of convolution feature map, and the calculation process of pooling feature map can be expressed as follows:(15)anl=fknl×1s2∑s×sanl-1+bnl,where s is the template size of downsampling and kml is the template weight. According to different sampling method, the pooling computing methods include maximum pooling, average pooling, and random pooling. Pooling is the aggregation statistics of the feature in a continuous region, and the regional feature is represented by the maximum value and average value and so on. In this paper, the maximum pooling is used.


## 2.4.3. Full Connection Layer

The full connection layer adopts the full connection mode, and the vector transformation of pooling feature map is processed. The two-dimensional feature matrix of this layer is stretched into one-dimension feature vector, which is convenient for the subsequent output layer calculation. The full connection layer is still equivalent to the convolution layer, and the difference is that the arrangement of the output is different, or the elements of pooling feature map are arranged directly into one-dimension vector.


## 2.4.4. Classifier Layer

The classifier layer mainly implements two functions: (1) In training stage of CNN, it provides the labels of the training data, and the training error is calculated. The weights are adjusted by gradient descent and reverse propagation, and the training error is reduced to minimum which completes the training process of CNN. (2) In testing stage of CNN, the output probability of testing samples is calculated, and the state of testing samples is determined according to the criterion of maximum probability. At present, the most commonly used classification method in classifier layer is logistic method and Softmax method. Because the degradation state recognition of planetary gear is a multipattern classification problem, in this paper, the Softmax method is used.


## 3. Experiment Introduction

The simulation experiment of the degradation state of planetary gear is carried out on the comprehensive simulation test bench for mechanical fault, and the basic structure of the comprehensive simulation test bench for mechanical fault is shown as Figure3.Figure 3
Simulation experiment of degradation states for planetary gear. (a) Controllable motor, (b) planetary gearbox, (c) vibration sensors, (d) spur gearbox, (e) load system, and (f) data acquisition instrument.In this experiment, the degradation states of two types of planetary gear faults are simulated, and they are degradation states of broken planetary gear and the degradation states of pitting planetary gear, respectively. Normal planetary gear and different degradation states are shown in Figure4. Based on the basic parameters of the planetary gearbox of test bench and the preliminary analysis of its vibration signal in the early stage, the feature frequency and its side frequency of planetary gearbox are in the frequency band with 20 Hz-640 Hz, and the highest and most prominent natural frequencies of planetary gearbox are in the frequency band with 2800 Hz-3200 Hz. According to the Nyquist sampling principle, the fault information contained in the frequency band of feature frequency should be focused, and the natural frequency of planetary gearbox also needs to be taken into account. In addition, considering that too high sampling frequency will increase the calculation amount and affect the calculation efficiency, and the sampling frequency of the experiment process is set as 6400 Hz. The basic parameters setting in the experiment are shown in Table 1.Table 1
Parameters setting of the experiment process.
Motor speed Sampling frequency Load Sample length 2700 r/min 6400 Hz 13.5 Nm 6400Figure 4
Normal planetary gear and different degradation states.


## 4. Experiment Analysis

The simulation experiment of degradation state of planetary gear is carried out, and the obtained vibration signals of normal planetary gear and different degradation states are shown in Figure5. Due to the nonlinear and nonstationary characteristics of the vibration signal of planetary gear and the tiny difference among the vibration signals of various degradation states, there is no obvious feature difference among vibration signals in time domain, and it is not possible to recognize the planetary gear state according to Figure 5.Figure 5
Vibration signals of normal planetary gear and different degradation states.Next, the Figure5 obtained vibration signals of planetary gear are processed to verify the validity of the proposed degradation state recognition method. First of all, the vibration signals are processed by SSD, in order to save the paper space, the SSD of the vibration signal of the degradation state of breakage level 2 is taken as an example. In the decomposition process, SSD can automatically realize the embedding dimension selection, and the frequency band of the decomposition result can be segmented actively by automatic component grouping, and the obtained SSDCs of the vibration signal of degradation state of breakage level 2 are shown in Figure 6.Figure 6
The obtained SSDCs of the vibration signal of degradation state of breakage level 2.It can be seen from Figure6 that the vibration signal with nonlinear and nonstationary characteristics of the degradation state of breakage level 2 is decomposed into 9 SSDCs by a series of iterations (in order to expression convenience, the final residual component is expressed as SSDC9), and SSDC1-SSDC9 are segmented from low frequency to high frequency. It is obvious that SSDC1-SSDC4 have a distinct periodicity, indicating that SSD can separate the periodic component hidden in the original vibration signal, and the high frequency information contained in the original vibration signal is separated into SSDC5-SSDC9. By SSD, each SSDC contains a large number of feature information that is benefit for recognizing degradation state of planetary gear. Next, the multiscale information dimension is used to extract and quantify the feature information from different time scales for each SSDC.In the calculation process of multiscale information dimension, the scale factorτ is set to 50. For each SSDC, the new time series under 50 time scales can be obtained, and then the information dimension features can be extracted for each new time series. It can realize the multidimensional feature expression from global signal information to detail signal information. The multiscale information dimensions of the SSDCs of normal planetary gear and different degradation states are shown as Figure 7.Figure 7
The multiscale information dimension of normal planetary gear and different degradation states. (a) Normal planetary gear, (b) breakage level 1, (c) breakage level 2, (d) breakage level 3, (e) breakage level 4, (f) pitting level 1, (g) pitting level 2, (h) pitting level 3, and (i) pitting level 4.
(a) (b) (c) (d) (e) (f) (g) (h)  (i)As shown in Figure7, for same SSDCs, the information dimensions extracted from different time scales have a certain difference, the values of information dimension show an upward trend from SSDC1-SSDC9, and, moreover, the values of information dimension also show a downward trend as the increasing of time scale. Because the vibration signal of normal planetary gear is relatively regular and simple, and when the fault occurs, the stiffness of planetary gear has a local nonlinear variation that results in vibration signal being nonlinear and nonstationary. The multiscale information dimension of normal planetary gear is smaller compared with that of the other malfunction degradation states. In addition, it can be seen from Figure 7 that the multiscale information dimensions of each SSDC are also different for different degradation states of planetary gear faults. The main difference is the information degradation state of planetary gear. For example, for breakage level 1, the information dimensions of SSDC7-SSDC9 are obviously smaller than those of other degradation states when the time scale is relatively large (time scale is greater 35). For breakage level 3, the information dimensions of SSDC1, SSDC4, and SSDC6 increase significantly when the time scale is near 17. For pitting level 1, the information dimensions of SSDC7 in time scale 5-10 are obviously different from that of other degradation states. Next, on the basis of the multiscale information dimension of each SSDC, CNN which can fully consider the internal relationship of the information dimension of adjacent SSDCs and adjacent scales are used to realize the effective recognition of the degradation state of planetary gear.The structure of CNN is built below, and the feature matrix constructed by the multiscale information dimension of each SSDC is defined as the input of CNN, so the feature matrix size of the input layer of CNN is 9 × 50. The hidden layer of CNN is composed of two convolution layers and two polling layers alternately. In convolution layer 1, the number of the convolution kernel is 6, and the size of convolution kernel is 3 × 3, the slip step is 1, and the activation function is selected as Sigmoid function. In polling layer 1, the size of polling area is 1 × 6, and the polling area is not overlapped. In convolution layer 2, the number of convolution kernel is 12, and the other parameters are the same as the first layer. In polling layer 2, the size of polling area is 1 × 3, and the polling area is not overlapped. In addition, the feature dimension of full connection layer is set as 120. Because 9 types of degradation states of planetary gear need to be recognized by CNN in this experiment, the various number of the output layer of CNN is set to 9, and the Softmax classifier is used. The basic parameters of CNN are shown in Table2.Table 2
Basic parameters of CNN.
Layers Input layer Convolution layer 1 Polling layer 2 Convolution layer 3 Polling layer 4 Full connection layer Classifier layer Feature matrix number 1 6 6 12 12 1 1 Feature matrix size 9 × 50 7 × 48 7 × 8 5 × 6 5 × 2 1 × 120 1 × 9 Kernel size  3 × 3 1 × 6 3 × 3 1 × 3The training of CNN is carried out, 200 training samples are selected for each kind of planetary gear state randomly, and there are 1800 training samples in all. The training samples are processed by SSD and multiscale information dimension, and the obtained feature matrixes are defined as the input of CNN for training. There are 9 neurons in the output layer of CNN, and the output vectors correspond to 9 types of degradation states of planetary gear. The training rate is set to 1, and the number of iterations is set to 100. The training process of CNN is shown in Figure8. It can be seen from Figure 8 that the mean square error of training samples tends to be stable after 65 iterations, and the training process of CNN is completed. Next, the recognition ability of CNN for degradation state of planetary gear is verified, and 100 testing samples are randomly selected for each type of planetary gear state, a total of 900 testing samples, and those are recognized by the trained CNN. The recognition result is shown in Figure 9. Meanwhile, in order to illustrate the advantage of CNN, other recognition methods, BP neural network [26] and SVM [27], are used to perform comparative analysis. And the training samples and testing samples which are the same as that applied to CNN method are used to train and test BP neural network and SVM, respectively. The recognition results of CNN, BP neural network, and SVM for testing samples are shown in Table 3.Table 3
The recognition results of CNN, BP neural network and SVM.
Normal Breakage level 1 Breakage level 2 Breakage level 3 Breakage level 4 Pitting level 1 Pitting level 2 Pitting level 3 Pitting level 4 Overall recognition CNN 100% 98% 97% 99% 94% 97% 95% 98% 97% 97.2% BP 90% 88% 91% 87% 83% 90% 86% 91% 89% 88.3% SVM 96% 93% 90% 97% 93% 94% 90% 94% 94% 93.4%Figure 8
The training process of CNN.Figure 9
The recognition result of CNN for testing samples.It can be seen from Figure9 and Table 3 that the proposed degradation state recognition method of planetary gear combining multiscale information dimensions of SSD and CNN can obtain better recognition results. The overall recognition rate of CNN is 97.2%, and for normal planetary gear, the recognition rate of CNN is 100%. The planetary gear state with the lowest recognition rate of CNN is the planetary gear with pitting level 1, and it can also reach 94%. However, when BP neural network and SVM are used to recognize the degradation state of planetary gear, it can be found that their recognition rates are significantly lower than that of CNN. The overall recognition rate of BP neural network is only 88.3%, and the planetary gear states with the highest recognition rate of BP neural network are breakage level 2 and Pitting level 3 and only reach 91%. With the application of SVM, the recognition rates of various planetary gear states have been improved to a certain extent, the overall recognition rate reaches 93.4%, but compared with the application of CNN, its recognition effect still has a certain gap. The experimental results show that the proposed method in this paper is a reliable, accurate, and effective method for the degradation state recognition of planetary gear.


## 5. Conclusions

A new degradation state recognition method of planetary gear based on multiscale information dimension of SSD and CNN is proposed in this paper. The SSD developed from SSA is suitable for processing the vibration signal with nonlinear and nonstationary characteristics generated by planetary gear, and the embedding dimension and component grouping segmentation can be determined. The original vibration signal can be converted into a series of SSDCs which are easy to be analyzed, and the feature information of planetary gear is contained in each SSDC. Aiming at a series of SSDCs, multiscale information dimension combining multiscale analysis and fractal information dimension which is a quantitative extraction method of feature information is studied, and the multidimensional feature expression for each SSDC from global signal information to detail signal information can be realized. The feature matrix composed of multiscale information dimension of each SSDC is defined as the input of CNN, and the training samples including different planetary gear states are used to train CNN, so that the trained CNN has the ability to recognize the degradation state of planetary gear effectively. The experimental results show that the proposed method is suitable for processing and analyzing the vibration signal of planetary gear, and the overall recognition rate of various planetary gear states is up to 97.2%. This method is an effective method for feature extraction and degradation state recognition of planetary gear.



---

*Source: 8716979-2019-03-10.xml*