A novel motor fault diagnosis method based on principal component analysis (PCA) with a discrete belief rule base (DBRB) system

The working environment of the motor is diverse and the vibration signal is typically non-stationary. Traditional filtering methods cannot be used to effectively analyze non-stationary signals which often also suffer from interference noise during the acquisition process, making fault feature extraction very difficult. Thus, it is important to eliminate noise from the signal. Wavelet threshold denoising [40] has been widely applied to signal denoising (figure 1).

Zoom In Zoom Out Reset image size

The threshold function is the key to the effect of wavelet threshold denoising. Traditional threshold functions include both soft and hard threshold functions, as follows:

$$\begin{equation}{\mathop w\limits^ \wedge} _{j,k} = \left\{ {\begin{array}{*{20}{l}} {sgn({w_{j,k}})(\left| {{w_{j,k}}} \right| - \tau ){} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} \left| {{w_{j,k}}} \right| \geqslant \tau } \\ {0{} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} \left| {{w_{j,k}}} \right| \geqslant \tau } \end{array}} \right.\end{equation} \tag{ 1 }$$

where ${w_{j,k}}$ is the wavelet coefficient; ${\mathop w\limits^ \wedge} _{j,k}$ is the wavelet coefficient after processing; τ is the threshold. sgn is a sign function, when ${w_{j,k}}$>0, sgn(${w_{j,k}}$) = 1; when ${w_{j,k}}$ = 0, sgn(${w_{j,k}}$) = 0; when ${w_{j,k}}$<0, sgn(${w_{j,k}}$) = −1.

The mathematical model of the hard threshold function is as follows:

$$\begin{equation}{\mathop w\limits^ \wedge} _{j,k} = \left\{ {\begin{array}{*{20}{l}} {{w_{j,k}}{} {} {} {} {} {} \left| {{w_{j,k}}} \right| \geqslant \tau } \\ {0{} {} {} {} {} {} {} {} {} {} {} {} {} {} {} \left| {{w_{j,k}}} \right| < \tau } \end{array}} \right..\end{equation} \tag{ 2 }$$

When statistical methods face multivariate problems, the complexity of the algorithm will be greatly improved. In many cases, if two variables are correlated, it means that there is overlapping information between them. To solve this problem, PCA recombines the original variables into a group of independent comprehensive variables. The central idea is to reduce the high-dimensional features to a few unrelated features and reflect the original information as much as possible [41–43]. Specific methods are as follows:

Assuming that the input array has m features and n sets of data, the output results in a two-dimensional topological structure with j neurons. The specific vector process is as follows:

• (a)  
  Input (n,m) eigenmatrix

  $$\begin{equation}X = ({x_{ij}}) = \left( {\begin{array}{*{20}{c}} {{x_{11}},}&{{x_{12}}, \cdots }&{{x_{1m}}} \\ {{x_{21}},}&{{x_{22}}, \cdots }&{{x_{2m}}} \\ \begin{gathered} \vdots \\ {x_{n1}}, \\ \end{gathered} &\begin{gathered} \\ {x_{n2,}} \cdots \\ \end{gathered} &\begin{gathered} \cdots \\ {x_{nm}} \\ \end{gathered} \end{array}} \right)\end{equation} \tag{ 3 }$$

  where sample feature ${x_{ij}}$ represents the jth feature of the ith sample in the sample.
• (b)  
  Correlation analysis

  $$\begin{align}\rho = \frac{\sum\limits_{i = 1}^n \left( x_{ij} - \frac{1}{n}\sum\limits_{i = 1}^n x_{ij}\right)\left(x_{ij + \tau } - \frac{1}{n}\sum\limits_{i = 1}^n x_{ij + \tau }\right) }{\sqrt {\sum\limits_{i = 1}^n {\left(x_{ij}^{} - \frac{1}{n}\sum\limits_{i = 1}^n x_{ij}\right)^2} } \sqrt {\sum\limits_{i = 1}^n {\left(x_{ij + \tau } - \frac{1}{n}\sum\limits_{i = 1}^n x_{ij + \tau }\right)^2} } }\end{align} \tag{ 4 }$$

  where $x_{ij + \tau }^{}$ represents the $j + \tau$th feature of the ith sample in the sample data, and $\rho$ is [−1,1], so as to represent the correlation of the two features. $\rho$ = −1 means that the features of the two vectors are negatively correlated, $\rho$ = 0 means that the two vectors are unrelated, and $\rho$ = 1 means that the two vectors are positively correlated.
• (c)  
  Standardized data

  $$\begin{equation}y_{ij}^{} = \frac{{nx_{ij}^{} - \sum\limits_{i = 0}^n {x_{ij}^{}} }}{{\sqrt {\frac{1}{{n - 1}}\sum\limits_{i = 1}^n {{{\left(nx_{ij} - \sum\limits_{i = 1}^n x_{ij} \right)}^2}} } }}\end{equation} \tag{ 5 }$$

  where $y_{ij}^{}$ is the standardized data sample, i is the sample number of input data, and j represents the characteristic dimension of input data.
• (d)  
  Normalized processing

  $$\begin{equation}Y = ({y_{ij}}) = \left( {\begin{array}{*{20}{c}} {{y_{11}},}&{{y_{12}}, \cdots ,}&{{y_{1m}}} \\ {{y_{21}},}&{{y_{22}}, \cdots ,}&{{y_{2m}}} \\ \begin{gathered} \vdots \\ {y_{n1}}, \\ \end{gathered} &\begin{gathered} \cdots \\ {y_{n2}}, \cdots ,\\ \end{gathered} &\begin{gathered} \vdots \\ {y_{nm}} \\ \end{gathered} \end{array}} \right)\end{equation} \tag{ 6 }$$

  where n and m are the corresponding dimensions of normalized data samples in equation (6), and $Y$ is the normalized feature matrix, aiming to map data results between (0,1).
• (e)  
  Solve the covariance matrix

  $$\begin{equation}S = \frac{1}{{n - 1}}Y{Y^T}\end{equation} \tag{ 7 }$$

  where, ${Y^T}$ is the transpose matrix of $Y$, S is the covariance matrix, and its eigenvalues and eigenvectors are solved to reduce its dimension. The non-negative eigenvalue is ${\lambda _1} \geqslant {\lambda _2} \geqslant \ldots \geqslant {\lambda _m} \geqslant 0$, and the orthogonal unit eigenvector is ${u_k}$.
• (f)  
  Principal component calculation

  $$\begin{equation}{Z_k} = {({u_k})^T}Y\end{equation} \tag{ 8 }$$

  $$\begin{equation}{v_k} = \frac{{{\lambda _k}}}{{\sum\limits_{k = 1}^m {{\lambda _k}} }}\end{equation} \tag{ 9 }$$

  where, ${Z_k}$ is the kth principal component (k ⩽ m), and ${v_k}$ is the variance contribution rate, which represents the proportion of the variance value of each principal component after dimension reduction to the total variance value. The larger the proportion, the more important the principal component is, and their sum is 1.

The core of BRB is ER, and the belief rule is expressed as follows:

$$\begin{align} & {\text{IF}} \left({x_1}is{A^k}_1\right){} {} {} {} {} {\text{and}}{} {} {} {} \left({x_2}is{A^k}_2\right)\ldots \left({x_{Tk}}is{A^k}_{Tk}\right)\nonumber\\ & {\text{THEN}}\left\{ {(Di{s_1},{\beta _{1k}}),(Di{s_2},{\beta _{2k}}),(Di{s_3},{\beta _{3k}})\ldots(Di{s_N},{\beta _{Nk}})} \right\},\nonumber\\ & \quad \left( {\sum\limits_{i = 1}^N {{\beta _{ik}} \leqslant 1} } \right)\nonumber\\ & {\text{with a rule weight}}{} {} {} {} {} {} {\theta _{\text{k}}}(k = 1,\ldots ,L)\nonumber\\ & {\text{and attribute weights}}{} {} {} {} {} {} {\delta _{\text{1}}},{\delta _2},\ldots ,{\delta _{Tk}} \end{align} \tag{ 10 }$$

where ${x_i}(i = 1,2,\ldots ,{T_k})$ represents the antecedent attribute, and $A_i^k$ represents the ith antecedent attribute's reference value, T_k is the antecedent attribute's number. $Di{s_n}(n = 1,2,\ldots ,N)$ is the consequent attribute's reference value, N is the reference value's number of the consequent attribute; ${\beta _{Nk}}$ is the reference value's confidence degree of consequent attribute; ${\theta _k}(k = 1,\ldots ,L)$ is the kth rule's weight; in the entire rule base, L is the total number of rules, ${\delta _i}$ is the antecedent attribute's relative importance.

The fusion method of BRB is based on ER rule, which mainly includes calculation of activation weight and output ER fusion. The overall algorithm flow is as follows:

• (a)  
  Determine the reference pointsThe linear piecewise function method is used to determine the degree of similarity between the input and the reference point. The formula is shown below

  $$\begin{align} {\alpha _{ij}} & = \frac{{u({a_{i(j + 1)}}) - x_i^{}}}{{u({a_{i(j + 1)}}) - u({a_{ij}})}}\nonumber\\[4pt] {\alpha _{i(j + 1)}} & = 1 - {\alpha _{ij}}{} {} {} {} {} {} \,\,\,if{} {} {} {} {} {} {} {} {} {} u({a_{ij}}) \leqslant x_i^{} \leqslant u({a_{i(j + 1)}})\nonumber\\[4pt] {\alpha _{ik}} & = 0{} {} {} {} {} {} {} {} {} {} {} for{} {} {} {} {} {} {} k = 1,\ldots ,{J_i},k \ne j,j + 1 . \end{align} \tag{ 11 }$$
• (b)  
  Calculate the each rule's activation weight

  $$\begin{equation}{w_k} = \frac{{{\theta _k} \times \prod _{i = 1}^{{T_k}}{{\left( {{\alpha _{ik}}} \right)}^{\overline {{\delta _i}} }}}}{{\sum\nolimits_{j = 1}^L {\left[ {{\theta _j} \times \prod\nolimits_{l = 1}^{{T_k}} {{{\left( {{\alpha _{lj}}} \right)}^{\overline {{\delta _l}} }}} } \,\right]} }}\end{equation} \tag{ 12 }$$

  where ${\bar \delta _{\text{i}}} = \frac{{{\delta _i}}}{\begin{gathered} \max \{ {\delta _i}\} \\ i = 1,\ldots {T_k} \\ \end{gathered} }$ represents the relative weight of the antecedent attribute, ${\alpha _{ik}}(i = 1,\ldots {T_k})$ represents the matching degree between ${x_i}$ and each reference value. Here,${\alpha _{ik}}$ ≧ 0, $\sum\nolimits_{i = 1}^{{T_k}} {{\alpha _{ik}}} \leqslant 1$, ${\alpha _k} = \prod \sum\nolimits_{i = 1}^{{T_k}} {({\alpha _{ik}}} {)^{{{\bar \delta }_i}}}$ is the joint matching degree. ${w_k}$ is the each rule's activation weight.
• (c)  
  ER fusionIn the previous step, the activation weight of each rule is calculated, and then the ER rule is used to combine multiple rules to generate the final conclusion, which refers to the confidence of the BRB model for each consequent attribute.

The confidence formula is as follows:

$$\begin{align} {\beta _j} & = \frac{{\mu \times \left[ {\prod\nolimits_{k = 1}^L \left({w_k}\beta _{jk} + 1 - {w_k}\sum\nolimits_{j = 1}^N \beta _{jk}\right) } \right] - \prod\nolimits_{k = 1}^L \left(1 - {w_k}\sum\nolimits_{j = 1}^N \beta _{jk}\right) }}{{1 - \mu \times \left[ {\prod\nolimits_{K = 1}^L {(1 - {w_k})} } \right]}}\nonumber\\ \mu & = \left[ \sum\nolimits_{j = 1}^N {} \prod\nolimits_{k = 1}^L \left({w_k}\beta _{jk} + 1 - {w_k}\sum\nolimits_{j = 1}^N \beta _{jk}\right) - (N - 1) \times \prod\nolimits_{k = 1}^L \left(1 - {w_k}\sum\nolimits_{j = 1}^N \beta _{jk}\right) \right]^{ - 1} .\end{align} \tag{ 13 }$$

Then the final fusion result is as follows:

$$\begin{equation}O(U) = \left\{ {(Di{s_j},{\beta _j}),{} {} {} {} j = 1,\ldots ,N} \right\}\end{equation} \tag{ 14 }$$

where $Di{s_j}$ is the output attribute of the jth consequent attribute, ${\beta _j}$ is the output confidence of the jth consequent attribute, k is the kth rule.

The establishment of the BRB needs to cover all the antecedent attributes and the each attribute's reference points. Therefore, in the BRB model, when T antecedent attributes exist and the number of reference points of the ith antecedent attribute is Ji, the BRB's rules can be represented as $\prod\limits_{i = 1}^T {{J_i}}$. For example, when the antecedent attribute's number is 5 and the number of reference points of each antecedent is 3, the rule numbers (RNs) of the BRB can be counted as 3⁵ = 243. Obviously, the RN in the BRB shows exponential growth. Therefore, a large number of antecedent attributes and reference points will lead to the 'combination explosion' problem, the most direct and effective solution is to reduce the number of antecedent attributes. To solve this problem, numerous studies have been carried out previously [44]. Reduced the number of antecedent attributes through grey target theory, multidimensional scaling, isomap and PCA. By evaluating the overall capability of armored vehicles, it was found that the BRB model based on PCA attribute reduction gave the best performance. However, when screening key attributes, the above method only considers the contribution degree of principal components and ignores the relative importance of each principal component. The amount of information of each principal component is different, so the importance of different principal components is different. The researchers used them as inputs for equally important traits, which distorted the original traits and rendered the diagnosis invalid.

Therefore, in this study it was proposed that the variance contribution rate of each principal component was taken as the relative weight of DBRB antecedent attributes, and the original characteristic information was restored. This not only makes up for the shortcomings of PCA, but also gives full play to the advantages of BRB, which can combine historical data and expert experience to deal with uncertainties and improve the operation speed and efficiency of the BRB model.

Wavelet soft threshold and wavelet hard threshold denoising are used to process the motor vibration signal respectively. The wavelet base is a db3 wavelet, the wavedec function in Matlab is used for multi-scale one-dimensional wavelet decomposition, and the wdencmp function is used for denoising. To compare the denoising effects of two different thresholds, Signal to Noise Ratio (SNR) and Mean Square Error (MSE) are used as the criteria to evaluate the denoising performance.

$$\begin{equation}{\text{SNR}} = 10\lg \frac{{\sum {f^2}(t)}}{{\sum {{(\mathop f\limits^ \wedge (t) - f(t))^2}}}}\end{equation} \tag{ 15 }$$

$$\begin{equation}{\text{MSE}} = \frac{1}{N}\sum {(\mathop f\limits^ \wedge (t) - f(t))^2}\end{equation} \tag{ 16 }$$

where $f(t)$ is the original signal, $\mathop f\limits^ \wedge (t)$ is the estimated signal after wavelet denoising, N is the length of the signal, t is the variable time, the sum operation along variable t. The larger the SNR and the smaller the MSE after denoising, the better the denoising effect.

Previous studies focused more on feature extraction of typical time domain features. In this study, typical time domain feature parameters and mathematical statistics feature parameters are used as dimension reduction and denoising feature parameters of PCA, and one-dimensional vibration signals are treated as a large number of values to carry out information mining using mathematical statistics to extract features of the data more effectively. The main feature parameters are shown in tables 1 and 2.

PCA is a useful statistical analysis method, which allows dimension reduction and denoising for the above 13 time domain mathematical statistical features. It reduces the above 13 features to three principal components with a high degree of combinatorial information.

Step 1: Determine the referential points

In order to avoid a large-scale rule base and a subsequently high calculation cost of the BRB, there should not be too many features corresponding to antecedent attributes. The referential points determined by expert experience are stochastic, and the division of referential points greatly affects the calculation speed and accuracy, so it is necessary to find the appropriate referential values. Here, a K-means algorithm [45] is used to determine the cluster centers, which are taken as the referential points, and the maximum and minimum values of features are taken as the reference maximum and minimum boundaries. The specific principle of the K-means algorithm is as follows:

• (a)  
  Cluster centers $\left\{ {{C_1},{C_2},{C_3}, \ldots ,{C_k}} \right\}$ is initialized, then the Euclidean distance from each object to the center of each cluster is calculated, its representation is as follows:

  $$\begin{equation}dis({X_i},{C_j}) = \sqrt {\sum\limits_{t = 1}^m {{{({X_{it}} - {C_{jt}})}^2}} } \end{equation} \tag{ 17 }$$

  where ${X_i}$ denotes the ith object, $1 \leqslant t \leqslant m$;${C_j}$ denotes the jth clustering center; ${X_{it}}$ denotes the tth attribute of the ith object;${C_{jt}}$ denotes the tth attribute of the jth cluster center; $1 < k \leqslant n$,$1 \leqslant j \leqslant k$, $1 \leqslant t \leqslant m$.
• (b)  
  The distances between each object and each cluster center are compared in turn, then the object is assigned to the class family of the closest cluster center to obtain K class groups, $\left\{ {{S_1},{S_2},{S_3}, \ldots ,{S_k}} \right\}$
• (c)  
  In each dimension of a class family, the average of all its objects is called the class center, and the formula is as follows:

  $$\begin{equation}{C_l} = \frac{{\sum\nolimits_{{X_i} \in {S_l}} {{X_i}} }}{{\left| {{S_l}} \right|}}\end{equation} \tag{ 18 }$$

  where ${C_l}$ denotes the lth cluster center,$\left| {{S_l}} \right|$ denotes the object's number in the lth family of classes, $1 \leqslant l \leqslant k$,$1 \leqslant i \leqslant \left| {{S_l}} \right|$.

Step 2: Establish discrete belief rule base (DBRB) model

After the reference values are determined, the OBRB calculates the match degree between the two reference values. To simplify the rule base, speed up diagnosis and reduce over-fitting, the DBRB is proposed. When the input is discrete, the interval of every two reference values is classified as one, and the similarity is directly assigned to 1, which produces a reference point.

Step 3: Calculate the matching degree of input data

The confidence degree of each category into which the feature falls is directly assigned as 1, and that of the others is assigned as 0. The corresponding confidence degree of the three categories are [1 0 0], [0 1 0], [0 0 1] respectively.

Step 4: Calculate the activation weights of input samples obtained online for corresponding rules.

Step 5: The analytic ER algorithm is used to fuse the activated rules.

Step 6: Optimize the established model

Through steps 1–5, the DBRB model is established. However, as the initial parameters are determined by experts, they have certain randomness and limitations, so optimal diagnosis results cannot be obtained. Therefore, we use historical data to train parameters and introduce an optimization method to optimize the DBRB model. Typically, parameters that can be optimized are output belief degree ${\beta _{jk}}$, initial rule weight ${\theta _k}$, and input weight ${\delta _T}$. To optimize the model, the objective function should be built first. The optimization system is shown in figure 3.

Zoom In Zoom Out Reset image size

P is the parameter to be optimized, which is the vector of the above three adjustable parameters:

$$\begin{equation}P = ({\beta _{ik}},{\theta _k},{\delta _j};{} {} {} {} {} {} {} i = 1,\ldots ,N;{} {} {} {} {} {} k = 1,\ldots ,L;{} {} {} {} {} {} {} j = 1,\ldots ,T).\end{equation} \tag{ 19 }$$

The number of optimization parameters p is (N *L +L +T). The fmincon function in MATLAB is selected to seek the optimal parameters' value in the specified interval. $\mathop {({\textrm O})}\limits^ \wedge$ is the form of observed belief degree for the real system; $({\textrm O})$ is the estimated belief degree set of system which is obtained through simulation. The minimum Euclidean distance (MED) of the observed belief degree $\mathop {({\textrm O})}\limits^ \wedge$ and the estimated belief degree $({\textrm O})$ are selected as the objective function of optimization. Where $\zeta (p)$ is the objective function to be optimized.

$$\begin{align} & \mathop {\min }\limits_p {} {} {} {} {} {} {} \left\{ {\zeta (p)} \right\}\nonumber\\ & \zeta (p)\, = \,\frac{1}{n} \times \sum\limits_{i = 1}^n {\sqrt {{{(\mathop {\textrm O}\limits^ \wedge - {\textrm O})}^{\text{2}}}} }\nonumber\\ & s.t.{} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} {} 0 \leqslant {\beta _{ik}} \leqslant 1 \nonumber\\ & \sum\limits_{i = 1}^3 {{\beta _{ik}}} = 1 \nonumber\\ & {\text{0}} \leqslant {\theta _k} \leqslant 1 \nonumber\\ & 0 \leqslant {\delta _j} \leqslant 1 . \end{align} \tag{ 20 }$$

To verify the ability of the proposed algorithm to diagnose faults, real fault data was collected. The relevant data samples of the vibration signals in this paper were collected from a permanent magnet motor TZ205XS70K01 provided by UNIlink Energy Innovation company, (figure 4). The motor parameters are shown in table 3.

Zoom In Zoom Out Reset image size

Three one-dimensional vibration acceleration signals were collected at a sampling frequency of 20 kHz. The three working conditions were: G1—bearing inner ring fault, G2—rotor eccentricity, G3—stator short circuit. A total of 300 000 data samples for these three conditions were extracted separately from the test bench.

Wavelet soft threshold denoising and wavelet hard threshold denoising methods are used to obtain an improved denoising signal. For intuitive display, the results of 6000 data points in a time window are shown in figure 5. The abscissa is the sample, and the ordinate is the vibration acceleration.

Zoom In Zoom Out Reset image size

Obviously, the vibration signal has different degrees of improvement after the wavelet threshold denoising. The amplitude fluctuation becomes smaller and the curve is smooth. To compare the denoising effect, SNR and RMSE were used as reference indexes. With the wavelet soft threshold algorithm, SNR = 12.5701 and RMES = 0.2182; With the wavelet hard threshold algorithm, SNR = 10.4368 and RMES = 0.2847.

The rated speed of the motor is 3000 r min⁻¹ and the sampling frequency is 20 kHz. Each revolution is treated as a cycle, and each cycle contains 400 data samples. In order to select an appropriate time length to characterize each time domain information and avoid the randomness of small-period signal fluctuations, 15 cycles were selected as a data sample. Every 6000 data points are taken as a time window in chronological order, there are 50 data samples in each working condition data, a total of 6000 × 50 data sets. Therefore, three working conditions constitute a total of 150 data features, thus forming a data matrix of (150,6000). Each working condition is named label 1, label 2, and label 3 respectively.

After the original one-dimensional data was segmented into time series, a time matrix of 150 rows and 6000 columns was formed, and features were extracted from each row of data. A total of seven typical time domain features and six digital statistical features were extracted.

The 13 extracted features are reduced to three principal components, and the feature visualization is shown in figure 6 below. After PCA dimension reduction, all features are separated, with only one feature of G2 and G3 slightly overlapping, which reduces the complexity of fault identification. The contribution rate of variance of the first component is about 82%, that of variance of the second component is about 12%, and that of variance of the third component is 6%. The complete data processing flow diagram is shown in figure 7:

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Before establishing the DBRB model, the precursor properties are set, ${\delta _{\text{1}}} = {v_1} = 0.82$, ${\delta _2} = {v_2} = 0.12$, ${\delta _3} = {v_3} = 0.06$. The variance contribution rate is calculated according to Formula 9. In order to simplify the rule base, the clustering quantity K is set as 2,1,1 respectively, so that the first feature clustering center is −9.014, 51.242; The clustering center of the second feature was 11.470; The clustering center of the third feature was 0.661, and the maximum and minimum of each feature is taken as the reference value boundary. The referential points are as follows:

$$\begin{equation*}\begin{array}{*{20}{l}} {{A_1} = \left\{ { - 93.323, - 9.014,51.242,98.85} \right\}} \\ {{A_2} = \left\{ { - 52.646,11.470,50.833} \right\}} \\ {{A_3} = \left\{ { - 4.151,0.661,5.776} \right\}} \end{array}.\end{equation*}$$

If traditional continuous BRB (OBRB) is used here, the initial rule base should include 36 (4 × 3 × 3) rules. In contrast, DBRB is introduced to reduce the referential points, where Ll means Low level, ml means Middle level, Hl means High level. So the discrete interval is divided as follows:

$$\begin{equation*}\begin{array}{*{20}{l}} {{A_1} = \left\{ {Ll = 1;Ml = 2;Hl = 3} \right\}} \\ {{A_2} = \left\{ {Ll = 1,Hl = 2} \right\}} \\ {{A_3} = \left\{ {Ll = 1,Hl = 2} \right\}} \end{array}.\end{equation*}$$

According to the interval form of 3, 2, and 2, the initial rule base of vector form is established as shown in table 4. It includes 12 rules (3 × 2 × 2), the size of rule base is reduced by 2/3. The initial rule weights are all set to 1.

Table 4. Initial DBRB.

No.	θk	X1 and X2 and X3	Confidence distribution structure
1	1	Ll and Ll and Ll	$\left\{ {\left( {Di{s_1},1} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},0} \right)} \right\}$
2	1	Ll and Ll and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},1} \right),\left( {Di{s_3},0} \right)} \right\}$
3	1	Ll and Hl and Ll	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},1} \right)} \right\}$
4	1	Ll and Hl and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},1} \right),\left( {Di{s_3},0} \right)} \right\}$
5	1	Ml and Ll and Ll	$\left\{ {\left( {Di{s_1},1} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},0} \right)} \right\}$
6	1	Ml and Ll and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},1} \right)} \right\}$
7	1	Ml and Hl and Ll	$\left\{ {\left( {Di{s_1},1} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},0} \right)} \right\}$
8	1	Ml and Hl and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},1} \right),\left( {Di{s_3},0} \right)} \right\}$
9	1	Hl and Ll and Ll	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},1} \right)} \right\}$
10	1	Hl and Ll and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},{\text{1}}} \right),\left( {Di{s_3},{\text{0}}} \right)} \right\}$
11	1	Hl and Hl and Ll	$\left\{ {\left( {Di{s_1},1} \right),\left( {Di{s_2},0} \right),\left( {Di{s_3},0} \right)} \right\}$
12	1	Hl and Hl and Hl	$\left\{ {\left( {Di{s_1},0} \right),\left( {Di{s_2},1} \right),\left( {Di{s_3},0} \right)} \right\}$

5.6.1. Fault diagnosis using the DBRB model.

In this case, three principal components are taken as three features, each fault has 300 000 vibration datasets, and 6000 datasets are taken as a time window. Then each working condition has 50 data segments, that is, each feature has 50 data segments, and each feature has a total of 150 samples for three working conditions. Here, 100 samples and 50 samples are randomly selected as the training set and the test set, respectively. After optimization of the DBRB model, the results of the training set are shown in figure 8(a) and the results of the test set are shown in figure 8(b).

Zoom In Zoom Out Reset image size

In the figure 8(a), the accuracy of the training set is 92%. The accuracy of the test set is up to 96%, and the MED is 0.7146. The optimized D-BRB structure is shown in table 5. Attribute weights before and after optimization are shown in table 6.

Table 5. Optimized DBRB structure.

No.	θk	X1 and X2 and X3	Confidence distribution structure
1	1	Ll and Ll and Ll	$\left\{ {\left( {Di{s_1},0.024} \right),\left( {Di{s_2},0.953} \right),\left( {Di{s_3},0.023} \right)} \right\}$
2	1	Ll and Ll and Hl	$\left\{ {\left( {Di{s_1},0.971} \right),\left( {Di{s_2},0.012} \right),\left( {Di{s_3},0.017} \right)} \right\}$
3	1	Ll and Hl and Ll	$\left\{ {\left( {Di{s_1},0.013} \right),\left( {Di{s_2},0.013} \right),\left( {Di{s_3},0.974} \right)} \right\}$
4	1	Ll and Hl and Hl	$\left\{ {\left( {Di{s_1},0.013} \right),\left( {Di{s_2},0.973} \right),\left( {Di{s_3},0.014} \right)} \right\}$
5	1	Ml and Ll and Ll	$\left\{ {\left( {Di{s_1},{\text{0}}{\text{.941}}} \right),\left( {Di{s_2},0.{\text{029}}} \right),\left( {Di{s_3},0.{\text{030}}} \right)} \right\}$
6	1	Ml and Ll and Hl	$\left\{ {\left( {Di{s_1},0.017} \right),\left( {Di{s_2},0.017} \right),\left( {Di{s_3},0.966} \right)} \right\}$
7	1	Ml and Hl and Ll	$\left\{ {\left( {Di{s_1},{\text{0}}{\text{.012}}} \right),\left( {Di{s_2},0.{\text{975}}} \right),\left( {Di{s_3},0.{\text{013}}} \right)} \right\}$
8	1	Ml and Hl and Hl	$\left\{ {\left( {Di{s_1},0.955} \right),\left( {Di{s_2},0.019} \right),\left( {Di{s_3},0.026} \right)} \right\}$
9	1	Hl and Ll and Ll	$\left\{ {\left( {Di{s_1},0.012} \right),\left( {Di{s_2},0.012} \right),\left( {Di{s_3},0.976} \right)} \right\}$
10	1	Hl and Ll and Hl	$\left\{ {\left( {Di{s_1},0.023} \right),\left( {Di{s_2},0.953} \right),\left( {Di{s_3},0.024} \right)} \right\}$
11	1	Hl and Hl and Ll	$\left\{ {\left( {Di{s_1},0.954} \right),\left( {Di{s_2},0.021} \right),\left( {Di{s_3},0.025} \right)} \right\}$
12	1	Hl and Hl and Hl	$\left\{ {\left( {Di{s_1},0.{\text{019}}} \right),\left( {Di{s_2},{\text{0}}{\text{.013}}} \right),\left( {Di{s_3},0.{\text{968}}} \right)} \right\}$

Table 6. Attribute weights comparison.

	Before optimization			After optimization
Weight	${\delta _{\text{1}}}$	${\delta _{\text{2}}}$	${\delta _{\text{3}}}$	${\delta _{\text{1}}}$	${\delta _{\text{2}}}$	${\delta _{\text{3}}}$
Value	0.82	0.12	0.06	0.82	0.26	0.13

Interestingly, the first antecedent weight${\delta _{\text{1}}}$ was the same before and after optimization, the first principal component information is preserved completely, which further illustrates the superiority of the proposed method in defining the weight factor.

Here, the variance contribution rate is set as the initial attribute weight in order to better express and restore the initial feature information, so as to better show the real operation state of the motor. In this paper, the initial attribute weights are optimized by DBRB model in order to find the parameters that can best fit the model, so as to get the optimal diagnosis results.

5.6.2. The influence of feature number on diagnosis result.

To explore the influence of feature numbers on diagnosis accuracy, mathematical statistical features were gradually added while traditional time domain features were retained. 12, 11, 10, 9, 8 and 7 features were selected respectively and inputted into DBRB model after PCA dimension reduction. The average value was taken after 20 experiments. The results of diagnosis accuracy and MED are shown in figures 9(a) and (b).

Zoom In Zoom Out Reset image size

As can be seen from the figures above, when the feature numbers are 10 and 11, although the diagnostic accuracy is the same, the MED increases. With a decrease in feature number, the diagnosis accuracy gradually decreases and the MED gradually increases. Therefore, it is necessary to select mathematical statistical features as the initial feature extraction parameters, which increases the amount of information containing data to deliver improved fault diagnosis.

5.6.3. DBRB compared with OBRB.

5.6.4. DBRB compared with other popular classification algorithms.

In a similar way, the DBRB model developed here was compared with several traditional mainstream classifiers, the three features after feature dimension reduction are taken as input attributes, 50, 40 and 30 samples of 150 data samples are selected as test sets, and the average value of 30 experiments is taken respectively. The Naive Bayes (NB), Bayes Net (BN), Random Forest (RF) and Support Vector Machine (SVM) classifiers were selected for comparison. Their diagnostic accuracy is summarized in figures 10(a)–(c). The upper and lower edges of the box are the upper quartile (Q25%) and the lower quartile (Q75%) of the accuracy, the upper and lower edges represent the maximum and minimum values of the data, the small rectangle is the average value of the accuracy, and the black diamond outside the box is the outlier.

Zoom In Zoom Out Reset image size

As can be seen from figure 10, the DBRB has the highest diagnostic accuracy compared to other algorithms, with no outliers and a more concentrated accuracy distribution, indicating that the DBRB model is more robust. To further verify the reliability of experimental results, 40 test samples and 30 test samples were used in the training set and test set. Again, similar conclusions can be drawn, confirming the superior accuracy and robustness of the presented model for fault diagnosis from inherently uncertain motor vibration signal data.