An effective multi-channel fault diagnosis approach for rotating machinery based on multivariate generalized refined composite multi-scale sample entropy

Fault diagnosis of critical rotating machinery components is necessary to ensure safe operation. However, the commonly used rotating machinery fault diagnosis methods are generally based on the single-channel signal processing method, which is not suitable for processing multi-channel signals. Thus, to extract features and carry out the intelligent diagnosis of multi-channel signals, a novel method for rotating machinery fault diagnosis is proposed. Firstly, a novel nonlinear dynamics technique named the multivariate generalized refined composite multi-scale sample entropy was presented and applied to extract fusion entropy features of multi-channel signals. Secondly, a practical manifold learning known as supervised isometric mapping was introduced to map the high-dimensional fusion entropy features in a low-dimensional space. In a third step, the Harris hawks optimization-based support vector machine was applied to carry out the intelligent fault recognition. Finally, aiming to verify the effectiveness of the proposed method and present its advantages, it was applied to analyze the rotating machinery system bearing and gear data. The experimental results have shown that the method at hand can accurately identify various faults in both the bearings and gears. Furthermore, in addition to being suitable for multi-channel signal fault diagnosis, it had higher recognition accuracy compared to other multi-channel or single-channel methods.


Introduction
Rotating machinery has been widely employed in various industrial fields, including petroleum, aerospace, metallurgy, energy, and chemistry. As critical rotating machinery components, bearings and gears have a significant role when it comes to supporting shafts, power conversion, and power transmission. However, due to complex operating conditions, their performance gradually degrades during a long-term operation and faults occur. Such behavior affects the performance of rotating machinery, their reliability, and service life, often causing inevitable economic losses. Thus, it is necessary to implement an effective intelligent fault diagnosis method for both bearings and gears [1][2][3].

Motivation
During the previous two decades, intelligent fault diagnosis methods based on vibration signals were gradually developed and applied in rotating machinery. Such methods were commonly used as they were proved to be effective [4]. They mostly included two steps: (1) the signal processing algorithms were firstly used to extract fault features and (2) various pattern recognition algorithms were applied to diagnose the faults based on the extracted features [5]. Commonly used intelligent fault diagnosis methods were often based on single-channel signals [6][7][8][9][10].
However, for rotating machinery, a single information source often contains limited fault information, meaning that it may result in ambiguous or even plain wrong diagnosis. Additionally, with the continuous improvement of rotating machinery safety requirements, this type of technology cannot meet the growing demand for data and information mining [11]. Thus, the fusion data processing methods based on multi-channel signals attracted the attention of scholars aiming to detect and diagnose the situation in detail. Such multi-channel signal processing methods are slowly becoming a standard part of intelligent fault diagnosis development [12,13]. For these technologies, signals obtained by the multiple sensors are processed as a whole, enabling them to fully extract the associated information, improve the information utilization rate, and finally obtain more accurate and reliable diagnosis results.
Moreover, the fusion features obtained using the multi-channel rotating machinery signals are regularly high-dimensional, nonlinear, and mixed with redundancy. Such properties increase the recognition time and reduce the diagnostic accuracy of subsequent classifiers [14]. When it comes to selecting lowdimensional and sensitive feature sets from these massive and nonlinear fusion features, it is imperative to consider intelligent fault diagnosis methods.

Entropy-based feature extraction approach
The rotating machinery vibration signals are usually nonlinear and non-stationary, which significantly increases the fault feature extraction complexity [14]. Fortunately, with the continuous development of nonlinear theory, entropy-based feature extraction algorithms which can overcome this problem were developed, as mentioned earlier. As such, they are widely applied for the feature extraction of various mechanical signals. The said technology mainly includes approximate entropy, sample entropy, permutation entropy, fuzzy entropy, and dynamic symbolic entropy [15][16][17][18][19].
Sample entropy is a frequently used method, which can effectively detect the minuscule time series changes. When the rotating machinery faults occur, the corresponding nonlinear dynamics complexity changes accordingly, requiring the application of sample entropy-based feature extraction algorithms [20][21][22][23]. For example, Wang combined the sample entropy method with ensemble empirical mode decomposition to identify the rolling bearing faults [20]. Gao et al. combined the multi-scale sample entropy with local mean decomposition to effectively extract the wind turbine fault features [21]. Furthermore, Dai et al. proposed a composite multi-scale sample entropy algorithm developed by combining the composite coarse-grained structure and the multiscale sample entropy. The algorithm was proved to effectively extract the bearing features [22]. Recently, Wang et al. have proposed the generalized refined composite multi-scale sample entropy algorithm and verified its performance by using bearing data [23].
The above-mentioned sample entropy-based algorithms performed very well in rotating machinery fault diagnosis; however, all the algorithms belong to the group of single-channel analysis methods. Since such methods are not suitable for processing multi-channel signals, it is vital to explore the state-of-the-art multichannel signal processing methods. Lu and Wang proposed the multivariate multi-scale sample entropy (MMSE) algorithm based on multi-scale sample entropy [24]. The algorithm was used to analyze any number of channel signals, providing reliable complexity analysis results. Nevertheless, when using MMSE to extract the fusion features of multi-channel rotating machinery signals, two shortcomings remained. (1) The new multivariate coarse-grained structure (MCS) of MMSE only considers the single MCS under the scale factor s, ignoring other equally important MCSs. (2) MMSE obtains a new MCS by sliding a window and averaging all the data within, which can cause the amplitude information of multichannel signals to be ignored [25,26].
However, a number of studies were carried out aiming to improve the MMSE. To solve the first shortcoming, Humeau-Heurtier proposed the multivariate refined composite multi-scale sample entropy (MRCMSE) algorithm. The algorithm considered multiple MCSs (i.e., multivariate composite coarsegrained structures) under s. Furthermore, the refined operation was used to avoid invalid entropy values [25]. Additionally, regarding the second deficiency, Yin [26].
Unfortunately, both the MRCMSE and the MGMSE only solve one MMSE shortcoming. Thus, a novel MGRCMSE algorithm was proposed in this paper, aiming to overcome both defects by combining MRCMSE and MGMSE. To prove its effectiveness, the algorithm was applied to extract the fusion feature of multi-channel vibration signals.

The proposed multi-channel fault diagnosis method
In this study, the authors have proposed an intelligent fault diagnosis method based on the fusion entropy manifold features and optimized support vector machine. The proposed method enables the effective mining of the fusion features from multi-channel rotating machinery signals, increasing the fault diagnosis accuracy, which is the main contribution of this paper. The method is divided into three phases: • The first phase-a novel nonlinear algorithm, named multivariate generalized refined composite multi-scale sample entropy (MGRCMSE), was developed to enable mining the fusion features from multi-channel signals. • The second phase-includes utilizing the supervised isometric mapping (SIM) [23] algorithm to reduce the fusion entropy feature dimensionality, enabling the authors to obtain the low-dimensional and sensitive fusion entropy manifold feature sets. • The third phase-is focused on fault identification and diagnosis; the Harris hawks optimizationbased support vector machine (HHO-SVM) [27] was applied to both identify and diagnose the faults.
The effectiveness of the proposed method was verified by analyzing the fault diagnosis experiments of bearings and gears mounted within the rotating machinery system. Moreover, compared to other single-channel or multi-channel feature extraction methods, the proposed feature extraction method extracted features from the original signals more effectively and with a higher fault diagnosis ability.

Contributions and structure of this paper
The main contributions of this paper can be stated as follows: • A novel nonlinear dynamics analysis algorithm (MGRCMSE) was proposed. • The proposed algorithm was used to carry out the multi-channel fault diagnosis of the rotating machinery. • The proposed algorithm was combined with the SIM to enable sensitive feature extraction.
Finally, the paper was structured as follows: Sect. 2 outlines the mathematical background behind the MGRCMSE. The novel intelligent multi-channel fault diagnosis method for rotating machinery is presented in Sect. 3. The bearings and gears were used as test specimens in the experimental verification (Sect. 4). Finally, conclusions are given in Sect. 5.
2 Proposed multivariate generalized refined composite multi-scale sample entropy approach The MMSE and MGMSE algorithms were introduced in this section. Additionally, based on these two algorithms, the MGRCMSE was proposed after implementing the multivariate generalized composite coarse-grained structure and refined operation. The flowchart of these three algorithms allowing the easier comparison is shown in Fig. 1.

Multivariate multi-scale sample entropy
For a multi-channel time series X K ¼ fx K;i g N i¼1 , K = 1, 2, …, p is the number of channels, and N is the number of sample points of each channel time series. The specific MMSE steps are as follows: 1. The MCSs Y s K ¼ fy s K;j g Ns j¼1 of multi-channel signals at s are constructed using Eq. (1), where Ns ¼ ½N=s: 2. Similarly, for each Y s G; K , the similarity probabilities P s G; m ðrÞ and P s G; mþ1 ðrÞ can be calculated by applying the MMSE steps (2) to (5). 3. Lastly, the MGMSE output values can be computed as:

Multivariate generalized refined composite multi-scale sample entropy
Based on MGMSE, the MGRCMSE algorithm was proposed in this paper, which was developed combining the multivariate composite coarse-grained structure and refined entropy operation. For the key MGRCMSE steps are: 1. The multivariate generalized composite coarsegrained structures Y s G; C; K ¼ fy s G; C; K;j g Ns j¼1 of multi-channel signals at s are found using Eq. (9), where C = 1, 2, …, s: 2. For each Y s G; C; K , similarity probabilities P s G; C; m ðrÞ and P s G; C; mþ1 ðrÞ are calculated according to the MMSE steps (2) to (5).

The final MGRCMSE values can be computed by
using refined entropy operation:

Simulation analysis
The MGRCMSE parameter settings are determined by carrying out noise experiments, allowing the authors to compare its performance and effectiveness to the MMSE, MRCMSE, and MGMSE.

MGRCMSE parameter selection
In the MGRCMSE, four parameters including the data length N, embedding dimension m, similarity tolerance r, and scale factor s must be set manually. Moreover, to study the influence these parameters have on the MGRCMSE, three-channel white noise (wn) and three-channel pink noise (pn) were utilized.

Influence of the scale factor s on MGRCMSE.
There is no standard for adjusting the parameter s; therefore, in this paper, the authors set its value to s = 25. According to Fig. 2, the fusion entropy values of the three-channel white noise are higher than those of the three-channel pink noise for the majority of N values. The reason for this phenomenon is that, compared to pink noise, white noise has stronger randomness and weaker self-similarity; therefore, its entropy values are higher. Furthermore, the threechannel white noise (or three-channel pink noise) entropy curves are similar to each other regardless of N, meaning that they have a limited effect on the MGRCMSE. Thus, the noise signal length was set to N = 3000.
3. The influence of the embedding dimension m on MGRCMSE. By selecting four different m values as m = {1, 2, 3, 4}, the MGRCMSE was used to analyze two above-mentioned noise signals. The resulting entropy curves obtained using N = 3000, r = 0.15sd, and s = 25 are shown in Fig. 3. Figure 3 shows that, primarily, both the white and pink noise entropy values increase with the increase in m. However, the overall trend practically remains the same. Secondly, when m is a small value (for example, m = 1), the white and pink noise entropy curves have similar values on large scales. On the other hand, when m is a large value (i.e., m = 3 or 4), the entropy curves fluctuate noticeably. Finally, when m = 2, the MGRCMSE entropy curves of both noise signals are relatively smooth, making them more distinguishable on the larger scales. Therefore, the value m = 2 was used in this paper. 4. The influence of the similarity tolerance r on MGRCMSE. Four different r values (0.1sd, 0.15sd, 0.2sd, and 0.25sd) were used to find the optimal one; the MGRCMSE was applied to analyze the above-mentioned noise signals. Figure 4 illustrates the entropy curves for N = 3000, m = 2, and s = 25. Figure 4 shows that for the same noise signal, the MGRCMSE entropy curve decreases with the increase in r. This can be explained by the fact that, when the value of r is low, the number of matching templates increases, increasing the entropy values. Moreover, when r has a small value (e.g., r = 0.1sd), the MGRCMSE entropy curves of both noise signals fluctuate slightly. On the other hand, the number of matching templates decreases when r is a large value, reducing the entropy values. When r is large (e.g., r = 0.25sd), the MGRCMSE entropy curves of both noise signals are relatively close in the large scales. When r is 0.15sd or 0.2sd, the MGRCMSE entropy curves of both signal groups are relatively smooth and the noise signals can be distinguished on large scales. Therefore, r = 0.15sd was used in this paper.
To conclude, within the paper, the MGRCMSE parameters are set to N = 3000, m = 2, r = 0.15sd, and s = 25.

Comparison of the various multi-channel signal analysis algorithms
In this section, the MGRCMSE algorithm was compared to the MMSE, MRCMSE, and MGMSE algorithms, aiming to assess its effectiveness. Four threechannel signals (as shown in Table 1) were analyzed using each of the algorithms, and the analysis results are shown in Fig. 5   Lastly, compared to the MMSE, MGMSE, and MRCMSE, the proposed MGRCMSE can not only obtain smooth and easily distinguishable entropy curves but also provide the analysis results for the actual composite situation. Thus, significant improvements were identified when using MGRCMSE when measuring the complexity of multi-channel signals.

The proposed intelligent rotating machinery fault diagnosis method
In this section, a novel intelligent fault diagnosis method for rotating machinery was proposed. Following the proposal, related algorithms such as SIM and HHO-SVM were described in detail.

The proposed fault diagnosis scheme
The authors aimed to provide a precise and reliable rotating machinery fault diagnosis process under different working conditions (i.e., a robust process). Thus, an intelligent multi-channel fault diagnosis method for the rotating machinery was proposed, primarily to enable the multi-channel signal fault diagnosis by using MGRCMSE, SIM, and HHO- SVM. The flowchart of the proposed method is shown in Fig. 6. The specific steps needed to carry out the process are summarized as follows: 1. Multi-channel signal acquisition. The multi-channel vibration acceleration signals of rotating machinery were collected through multiple sensors and in various working states. Then, the collected samples are randomly divided into training samples (10%) and testing samples (90%). 2. Fusion entropy feature extraction. The proposed MGRCMSE method was used to extract fusion entropy features from multi-channel signals, allowing the extraction of a nonlinear and highdimensional fusion feature set. 3. Dimensionality reduction. The SIM algorithm was introduced to map the MGRCMSE set to the low-dimensional space, obtaining a low-dimensional and sensitive feature set. 4. Pattern recognition. The training feature set and testing feature set were normalized to the interval [0, 1], respectively. The training set was used as input to the HHO-SVM, allowing it to search for the optimal parameters (c best , g best ), allowing the authors to obtain the optimized SVM model. Finally, the testing set was input to the trained SVM classifier to diagnose faults.

Supervised Isomap manifold learning method
In the next step, the SIM algorithm was utilized to reduce the MGRCMSE feature set dimensionality, aiming to construct the sensitive, low-dimensional, and high-distinguish feature set. For a given data set the key SIM procedures are as follows: where dðu i ; u j Þ is the Euclidean distance between u i and u j , Lu i is the u i category, r is the average value of dðu i ; u j Þ, and g is the weighted coefficient. 2. The neighborhood graph is constructed using the k-nearest neighbor method. If two sample points (e.g., u i and u j ) are adjacent to each other, there is a side connection with a side length d s ðu i ; u j Þ; otherwise, there is no side connection. 3. The shortest path between any two samples can be calculated using the Dijkstra algorithm; the path is approximately regarded as the geodesic distance between these samples. 4. The geodesic distance matrix is mapped to the low-dimensional space by using the multi-dimensional scaling algorithm [28], allowing us to obtain the final dimensionality reduction result.

Harris hawks optimization-based support vector machine
The SVM, based on statistical learning theory, is the most commonly used classifier. It not can only deal with small samples and nonlinear classification problems but also has good generalization ability. However, its performance is closely related to two important parameters-c and g. Thus, to improve the SVM performance, this section introduces the HHO-SVM classifier obtained after using the novel HHO algorithm to find the optimal values of SVM parameters. The HHO-SVM is a novel classifier (flowchart shown in Fig. 7) utilizing the Harris hawks optimization algorithm (for more details, see [29]) to optimize the key SVM parameters. The specific steps are shown next.
1. Preprocess the data Equation (12) is applied to normalize the input data set between [0, 1].
where, r 1 and r 2 are random numbers between (0, 1). 6. Update the Harris hawk locations. The position of each Harris hawk can be updated by using Eqs. (14) to (18). If the fitness value of the new position is lower than that of the previous position, the previous position is replaced with the new one; otherwise, the previous position is retained.
1. Survey phase (i.e., jEj ! 1). The position of each Harris hawks is updated by: where HðtÞ and Hðt þ 1Þ are the Harris hawk positions at t and t ? 1 iteration, respectively. H rabbit (t) is the prey location, H rand (t) is the random Harris hawk position in the current iteration, and H m ðtÞ ¼ ½ P N p i¼1 H i ðtÞ=N p is the average position of all Harris hawks positions in the current iteration, while r 3 , r 4 , r 5 , r 6 , and q are the random numbers between 0 and 1. 2. Development phase (i.e., 0 jEj\1). The variable r 7 is defined as the probability of the prey escaping and is a random number between 0 and 1. For example, r 7 \0:5 indicates that the prey has succeeded in escaping while r 7 ! 0:5 indicates that it has failed. This phase can be further into four possible situations: • Soft besiege (i.e., r 7 ! 0:5; jEj ! 0:5)-the position of each Harris hawk is updated by: Hðt þ 1Þ ¼ DHðtÞ À EjJH rabbit ðtÞ À HðtÞj DHðtÞ ¼ H rabbit ðtÞ À HðtÞ • Hard besiege (i.e., r 7 ! 0:5; 0 jEj\0:5)-the position of each Harris hawk is updated by: • Soft besiege with progressive rapid dives, (i.e., r 7 \0:5; jEj ! 0:5)-the position of each Harris hawk is updated by: where D = 2 is the positional dimension, S is a random vector of size 1 Â D, LFðxÞ is the Levy flight function, and l , c are the random number between 0 and 1. Finally, b ¼ 1:5 and C(Á) is a Gamma function. • Hard besiege with progressive rapid dives, (i.e., r 7 \0:5; jEj\0:5)-the position of each Harris hawk is updated using: Fig. 7 HHO-SVM flowchart 7. Determine whether the stop condition is met. When the simulation reaches the previously selected maximum number of iterations, the cycle is terminated. The final prey position (c best , g best ) is used as the output. Otherwise, step (2) is repeated until the stop criterion is met. 8. Construct the optimal SVM training model. The final prey position is implemented in the SVM prediction model. 9. Identify the testing set categories. The testing samples are input into the trained SVM prediction model for classification with the recognition results as outputs.

Experimental verification
The experiments were carried out using two case studies as examples: the first being concerned with bearings and the second with gears. The aim was to evaluate the performance and verify the effectiveness of the proposed fault diagnosis method. Moreover, several comparative experiments were carried out to evaluate the proposed method and verify its advantages. These include both the comparison with other multi-channel feature extraction methods and the comparison with several single-channel feature extraction methods.

Experimental rig
The experimental data were collected from the drivetrain dynamics simulator (Fig. 8). The experimental platform consists of the variable speed drive motor, programmable motor controller, planetary gearbox, parallel shaft gearbox, brake, and programmable magnetic controller. The planetary gearbox parameters are given in Table 2.

Rolling bearing data set
Four bearing operation states were considered in this experiment: normal, inner ring fault, outer ring fault, and ball fault, with each being outlined in Fig. 9. The working experimental rig conditions were selected as follows: condition 1-the driving speed of 25 Hz and no load; and condition 2-the driving speed of 35 Hz and the 2 V load. The acceleration sensors were used to collect the vibration signals in the X, Y, and Z directions, each obtained with the sampling frequency of 3000 Hz (please see Table 3). Finally, the threechannel time-domain bearing waveforms under various operation states are shown in Figs. 10 and 11.

Experimental results and analysis on condition 1
According to the proposed fault diagnosis process, the MGRCMSE method was applied to extract the fusion features from three-channel bearing signals operating under condition 1. The results for 800 bearing samples analyzed under four different states are plotted in Fig. 12. Several phenomena are found in Fig. 12. Firstly, the mean entropy values of the normal bearing signals are higher when compared to the fault bearing signals, regardless of the scale. Such behavior is caused by the self-similarity of normal signals, which is weaker than that of fault signals. Therefore, the normal state entropy values are higher than their fault state counterparts. Hence, the MGRCMSE algorithm can effectively monitor the occurrence of bearing failure. Secondly, the entropy values of the four operation states at most scales are ordered as follows E NOR-[ E BF [ E IRF [ E ORF , which is consistent with the realistic case. Both phenomena confirm that the MGRCMSE algorithm can identify the faults.
However, Fig. 12 shows that the IRF and BF entropy mean curves are rather similar. Moreover, the extracted feature set is high-dimensional and redundant; if directly input to the classifier as such, both the recognition time and error will be increased. For this reason, according to the proposed fault diagnosis method, a supervised manifold learning algorithm (SIM) was adopted to generate the corresponding lowdimensional set. The results are shown in Fig. 13. It should be pointed out that based on several experiments, the optimal SIM parameters are d = 3, k = 92, and g = 0.4.
In the 3D visualization results (see Fig. 13), it was shown that SIM could entirely separate the four types of samples. Additionally, it is evident that samples are grouped based on their state. This proves that the SIM algorithm can effectively obtain a low-dimensional and sensitive set, which is easily distinguishable.
When classifying the set, the above-presented fusion entropy manifold (i.e., MGRCMSE ? SIM) feature set was used as input for the HHO-SVM classifier. The HHO algorithm was first used to find the optimal SVM parameters (i.e., c best = 79.86, g best-= 7.62), allowing us to establish the optimal SVM prediction model. The training stage fitness curve is shown in Fig. 14a. Testing samples were input to the previously trained SVM classifier for pattern recognition; output results and the confusion matrix are shown in Fig. 14b and c, respectively.
As shown in Fig. 14, 720 testing samples were correctly identified with an average recognition rate of 100%. More importantly, it was proved that the proposed fault diagnosis method could effectively and accurately diagnose bearings in different states.

Experimental results and analysis on condition 2
To investigate the effect of operating conditions on the proposed fault diagnosis method, it was used to analyze the bearing data under operating condition 2.
The MGRCMSE method was first applied to extract the fusion features from three-channel bearing signals.
The results for a total of 800 analyzed bearing samples under four different states are plotted in Fig. 15. Figure 15 shows that when the bearing conditions change, the normal state entropy curve is higher than that of the fault states (at larger scale numbers). Moreover, the entropy values of four analyzed fault states at most scales are ordered as This experiment further demonstrates the advantages of the proposed fault feature extraction method. Subsequently, the SIM algorithm was used to reduce the dimensionality of the collected MGRCMSE feature set, as shown in Fig. 16. It should be pointed out that the optimal SIM parameters were d = 3, k = 200, and g = 0.35. As shown in Fig. 16, the SIM method can be used to effectively separate four bearing states and obtain both low-dimensional and sensitive fault features.
Finally, the HHO-SVM classifier was used to identify faults of the obtained MGRCMSE ? SIM set. The HHO algorithm was applied to find the optimal SVM parameters (i.e., c best = 33.44, g best-= 2.15) necessary for building the optimal SVM prediction model. The training stage fitness curve is shown in Fig. 17a. Testing samples were then inputted to the previously trained SVM classifier for pattern recognition; output results and testing sample confusion matrices are shown in Figs. 17b and c, respectively.
It can be concluded that 720 testing samples were correctly identified with an average recognition rate of 100%. More importantly, it was proved that the proposed fault diagnosis method could, effectively and accurately, diagnose bearings states in various working conditions.

Gear data set
In this experiment, four gear operation states were established: normal, tooth breakage fault, missing tooth fault, and tooth root fault (Fig. 18). The experimental working conditions were as follows: the driving speed of 25 Hz and no load. Three acceleration sensors were used to collect the vibration signals in X, Y, and Z directions under different operation modes, obtained at the sampling frequency of 3000 Hz. The detailed description is shown in Table 4, while the three-channel gear waveforms under each operation mode are plotted in Fig. 19.

Experimental results and analysis
According to the proposed fault diagnosis method, the MGRCMSE algorithm was utilized to extract the fusion features from three-channel gear signals, as shown in Fig. 20.
The following observations are made based on Fig. 20: (1) The mean entropy values of normal gear signals were higher than that of faulty gear signals at most of the scales. The main reason for such behavior was the occurrence of vibration signals in the normal state, which have shown strong randomness, while the vibration signals in faulty states had a prolonged impact. Thus, the entropy values of fault states were smaller than that of the normal state, indicating that the MGRCMSE can monitor the occurrence of gear failure. (2) The mean entropy curves of the gear signals under four states can be separated from each other, showing that the proposed method can effectively distinguish between the fault states.
The SIM was used to reduce the extracted MGRCMSE feature set dimensionality. (The result is plotted in Fig. 21.) It should be pointed out that the optimal SIM parameters were obtained experimentally and are as follows: d = 3, k = 193, and g = 0.38.
The 3D visualization results (shown in Fig. 21) show that SIM thoroughly separates the four-state samples, completely avoiding sample overlap. Thus, the SIM effectiveness in the dimensionality reduction of the gear feature set was considered as verified.
Furthermore, the obtained fusion entropy manifold (i.e., MGRCMSE ? SIM) feature set was used as input for the classification using HHO-SVM. The HHO algorithm was firstly employed to find the optimal SVM parameters (c best = 56.99, g best-= 24.84), allowing us to establish the corresponding optimal SVM prediction model. The training stage fitness curve is shown in Fig. 22a, while the testing sample output results and the confusion matrix are shown in Fig. 22b, c.  Figure 22 shows that the HHO-SVM correctly identified all the 720 testing samples. Thus, the average recognition accuracy of the proposed method for the gear data set was 100%, proving that it can both effectively and accurately diagnose various gear operation states.

Performance comparison
The advantages of the proposed MGRCMSE method were then further explored using the bearing fault diagnosis as an example; two comparative studies were carried out.
1. Comparison with other multi-channel feature extraction algorithm. In this section, the MGRCMSE algorithm performance was compared to the performance of the MGMSE algorithm; the results are shown in Fig. 23. The MGMSE algorithm parameters in this study are listed in Table 5.
According to Fig. 23a-d, the mean MGMSE and MGRCMSE entropy curves for each bearing state are relatively similar. However, it is noted that the mean MGMSE entropy curves have slight fluctuations in large scales, especially for bearing without faults. Next, Fig. 23e-h shows that the entropy values obtained using MGRCMSE have low standard deviations on most scales. Additionally, several random entropy curves were drawn for each bearing state signal, using both algorithms, as shown in Fig. 23i-l. When compared with the MGRCMSE, the MGMSE entropy curves have larger fluctuations. These phenomena prove that the proposed MGRCMSE algorithm can obtain more stable entropy results when compared to the MGMSE.
Furthermore, the feature sets obtained using the above-mentioned algorithms were input into the HHO-SVM for recognition (for results, please see Fig. 24).
As illustrated in Fig. 24  These single-channel feature extraction methods were used to extract features from the X, Y, and Z channel signals. The extracted results were then submitted to the HHO-SVM classifier for fault identification and the results are shown in Fig. 25. The algorithm parameters used to carry out the study are given in Table 6.
As shown in Fig. 25, the diagnosis accuracy of the same single-channel feature extraction algorithms (i.e., GMSE, GRCMSE, GCMFE, and GCMPE) varied significantly depending on the channel signals. For example, the GRCMSE algorithm had the highest average recognition accuracy for channel X signals (i.e., 81.25%), while it was the lowest for channel Y signals (i.e., 78.19%). The total difference was 3.06%. Moreover, the proposed MGRCMSE method had the Fig. 23 Mean entropy curves, standard deviation values, and random group entropy curves, using two feature extraction algorithms. a NOR mean curves, b IRF mean curves c ORF mean curves, d BF mean curves, e NOR standard deviation, f IRF standard deviation, g ORF standard deviation, h BF standard deviation, i NOR entropy curves, j IRF entropy curves, k ORF entropy curves, and l BF entropy curves

Comparison with similar fault diagnosis methods.
To validate the potential application of the proposed fault diagnosis method, it was compared to other previously published state-of-the-art methods [33][34][35][36][37][38], as shown in Table 7. The

Conclusions
In this paper, a novel nonlinear dynamics technique named the multivariate generalized refined composite multi-scale sample entropy was proposed. The aim was to enhance the complexity of measuring the performance of multi-channel signals. The proposed method was then applied to extract the features of multi-channel rotating machinery signals. Both the simulation analysis and rotating machinery fault diagnosis experiments confirmed the effectiveness of the proposed method, along with its advantages. The following conclusions were made: 1. The proposed MGRCMSE algorithm can effectively extract the fusion entropy features for both bearings and gears. 2. Compared to the existing multi-channel feature extraction algorithm (i.e., MGMSE) and singlechannel feature extraction algorithms (i.e., GMSE, GRCMSE, GCMFE, GCMPE), the MGRCMSE displayed the highest fault diagnostic accuracy. 3. The proposed intelligent multi-channel fault diagnosis method, based on MGRCMSE, can accurately diagnose both the bearing and gear faults.  [23] GRCMSE Single-channel method m = 2, r = 0.15sd, N = 3000, and s = 25 Yu et al. [31] GCMFE Single-channel method m = 2, r = 0.15sd, n = 2, N = 3000, and s = 25 Zheng et al. [32] GCMPE Single-channel method m = 6, k = 1, N = 3000, and s = 25 Wang et al. (proposed) MGRCMSE Multi-channel method m = 2, r = 0.15sd, N = 3000, and s = 25 When the rotating machinery has faults, the corresponding vibration signals will show regular transients. The proposed MGRCMSE method was used for rotating machinery fault diagnosis by exploring the complexity and regularity of multi-channel vibration signals. However, the MGRCMSE method parameters rely on manual or empirical settings. Thus, in the future, the authors aim to focus on mitigating this problem, and further apply it to diagnose similar rotating machinery faults.