A machine learning-based rolling bearing fault diagnosis method

Abstract

The application discloses a kind of based on machine learning's rolling bearing fault diagnosis method, belong to mechanical equipment fault diagnosis technical field, comprising the following steps: step S1: VMD parameter optimization;Step S2: effective IMF component screening;Step S3: feature extraction and dimension reduction;Step S4: fault diagnosis.The application utilizes variational mode decomposition to reduce the noise in vibration signal, for the problem of variational mode decomposition parameter selection, adaptive grey wolf algorithm is used to optimize;The data after noise reduction is extracted using permutation entropy, and the dimension of fault feature is reduced using nonlinear principal component analysis algorithm, effectively reducing the influence of noise in vibration signal on fault diagnosis accuracy, which can significantly improve the accuracy of rolling bearing fault diagnosis.

Description

Technical Field

[0001] This invention relates to the field of mechanical equipment fault diagnosis technology, and specifically to a rolling bearing fault diagnosis method based on machine learning. Background Technology

[0002] Rolling bearings, as key support components in rotating machinery transmission systems, maintain the normal operating position and rotational accuracy of shafts and are widely used in many important engineering fields. However, they are also vulnerable components in rotating machinery; rolling bearing failures account for approximately 30% of rotating machinery failures, making them widely recognized as one of the main causes of rotating machinery malfunctions. Rolling bearing failures significantly impact equipment operating performance and safety; therefore, accurate fault diagnosis of rolling bearings has important practical engineering significance.

[0003] Currently, bearing fault diagnosis methods are mainly divided into two categories: traditional methods and machine learning methods. Traditional fault diagnosis methods are mostly based on signal analysis and kinematic analysis, and their result judgment heavily relies on expert experience, which limits their applicability. Machine learning fault diagnosis methods, on the other hand, utilize machine learning network models to automatically and intelligently learn the features of input data and classify data types based on extracted high-dimensional abstract features, thereby determining the fault type. Among these, machine learning algorithms such as Support Vector Machines, Random Forests, and Neural Networks have been applied in the field of rolling bearing fault diagnosis due to their classification capabilities. However, due to the complexity of the working environment, the original vibration signals of bearings often contain a large amount of noise, which increases the difficulty of signal feature extraction and thus affects the accuracy of bearing fault diagnosis. Therefore, this invention proposes a machine learning-based rolling bearing fault diagnosis method. Summary of the Invention

[0004] The technical problem to be solved by this invention is: how to reduce noise and extract features from bearing vibration signals, and then perform accurate fault diagnosis, and provides a rolling bearing fault diagnosis method based on machine learning.

[0005] The present invention solves the above-mentioned technical problems through the following technical solution, and the present invention includes the following steps:

[0006] Step S1: VMD parameter optimization

[0007] Vibration signals under different states are input into the variational mode decomposition algorithm, and the adaptive gray wolf algorithm is used to optimize the variational mode decomposition parameters, which include the number of modes K and the penalty factor α.

[0008] Step S2: Screening of Effective IMF Components

[0009] The vibration signal is decomposed using the optimized variational mode decomposition parameters to obtain multiple IMF components. The correlation coefficient and permutation entropy of each IMF component are calculated, and the effective IMF components corresponding to the original vibration signal are selected by double thresholding.

[0010] Step S3: Feature Extraction and Dimensionality Reduction

[0011] The effective IMF components are input into the kernelized nonlinear principal component analysis algorithm. The principal components whose cumulative contribution rate exceeds the set threshold are selected, and the fault feature set of rolling bearing vibration data is constructed by combining time domain and frequency domain features.

[0012] Step S4: Fault Diagnosis

[0013] The particle swarm optimization support vector machine model is trained using the fault feature set, and the trained model is used for fault diagnosis of rolling bearings.

[0014] Furthermore, in step S1, the specific processing procedure of the adaptive gray wolf algorithm is as follows:

[0015] Step S11: Generate the initial gray wolf position using Latin hypercube sampling;

[0016] Step S12: Calculate the initial fitness and determine the leader based on the fitness, where the leader is α, β, or δ wolf;

[0017] Step S13: Enter the iteration loop;

[0018] Step S14: Output the optimal solution.

[0019] Furthermore, in step S11, the specific processing procedure is as follows:

[0020] Step S111: Define the parameter range matrix

[0021] Construct a parameter range matrix of dimension d×2:

[0022]

[0023] Where d represents the dimension of the optimization problem. Let K be the minimum and maximum values ​​of the k-th dimension parameter. There are two parameters to be optimized: the number of modes K and the penalty factor α. The value of d is 2.

[0024] Step S112: Generate LHS sampling matrix

[0025] For each dimension k, generate N equally spaced intervals:

[0026]

[0027] Where i = 1, 2, ..., N are interval indices, corresponding to N individual gray wolves;

[0028] Randomly sample a point within each interval to form a sampling sequence of dimension k:

[0029]

[0030] Perform the above operation on all dimensions to generate d independent sampling sequences;

[0031] Step S113: Random permutation and combination

[0032] Sampling sequence S for each dimension k Randomly arrange the sampled sequences of each dimension and combine them to generate the initial population matrix:

[0033]

[0034] Here, rand(i) represents the random permutation index of the i-th row.

[0035] Furthermore, in step S12, the fitness calculation formula is as follows:

[0036] Fitness=ω1·MSE+ω2·(-SNR)

[0037] Where ω1 and ω2 are weighting coefficients; MSE is the mean square error between the denoised signal and the original vibration signal; and SNR is the signal-to-noise ratio of the denoised signal.

[0038] Furthermore, in step S13, the specific processing procedure is as follows:

[0039] Step S131: Dynamically adjust the convergence factor τ

[0040] Update τ using a nonlinear decay strategy:

[0041]

[0042] Where t is the iteration number; τ initial T is the initial convergence factor; max γ is the maximum number of iterations; γ is the decay coefficient.

[0043] Step S132: Location Update

[0044] Calculate the distance to the prey:

[0045]

[0046] Updated Grey Wolf location:

[0047]

[0048] in, Let be the distance vector between the individual gray wolf and its prey; This represents the location of the prey, i.e., the current optimal solution. Let be the current position of the gray wolf individual at the t-th iteration, i.e., the candidate solution; The coefficient vector is calculated using the following formula:

[0049]

[0050]

[0051] in, [0,1] d Random vectors between;

[0052] Step S133: Perform adaptive mutation processing, as shown in the following formula:

[0053]

[0054] Where σ is the variable asynchronous length; [0,1] d Random vectors between This represents the latest position of the prey after adaptive mutation processing, i.e., the current solution.

[0055] Step S134: Recalculate fitness and update α, β, and δ wolves;

[0056] Step S135: Check the termination condition. If the termination condition is met, exit the loop. The termination condition is reaching the maximum number of iterations or the optimal fitness is less than the threshold.

[0057] Furthermore, in step S2, a dual-threshold screening is set, that is, a correlation coefficient threshold ρ is set simultaneously. threshold and permutation entropy threshold PE threshold Retain those that simultaneously satisfy ρ i ≥ρ threshold And PE i ≤PE threshold The IMF component, i represents the i-th IMF component, i.e. the effective IMF component obtained by filtering.

[0058] Furthermore, the formula for calculating the correlation coefficient threshold is:

[0059] ρ threshold =μ ρ +0.5σ ρ

[0060] The formula for calculating the permutation entropy threshold is:

[0061] PE threshold =μ PE -1.2σ PE

[0062] Where, μ ρ μ PE σ represents the correlation coefficient and the mean of the permutation entropy. ρ σ PE For the correlation coefficient and the standard deviation of the permutation entropy, respectively.

[0063] Furthermore, in step S2, the correlation coefficient is calculated using the following formula:

[0064]

[0065] Where x(t) is the original vibration signal with a length of N; IMF i (t) represents the i-th IMF component obtained from VMD decomposition; The mean of the original vibration signal and the i-th IMF component;

[0066] The formula for calculating permutation entropy is:

[0067]

[0068] Where m is the embedding dimension; p j It represents the probability of the j-th permutation pattern.

[0069] Furthermore, in step S3, the specific processing procedure of the kernelized nonlinear principal component analysis algorithm is as follows:

[0070] Step S31: Gaussian kernel calculation

[0071] The Gaussian kernel K(x) is calculated using the following formula. i x j ):

[0072]

[0073] Where ψ is the Gaussian kernel bandwidth parameter; ||x i -x j || 2 For sample x i With x j The square of the Euclidean distance;

[0074] Step S32: Approximate accelerated kernel matrix calculation

[0075] Step S321: Random sampling

[0076] From N samples, m base samples are randomly selected, where m << N, and m is a random number of samples. Approximate base sample size;

[0077] Step S322: Calculate the submatrix

[0078] Calculate the kernel matrix between the base samples:

[0079] W ij =K(x) i x j )

[0080] Where W∈R m×m i, j = 1, 2, ... m;

[0081] Calculate the kernel matrix between the base sample and other samples:

[0082] C il =K(x) i x m+l )

[0083] Where C∈R m×(N-m) i is the base sample index; l is the non-base sample index;

[0084] Step S323: Approximate kernel matrix

[0085] The kernel matrix is ​​approximated using the following expression:

[0086]

[0087] in, This is an approximation matrix of the original kernel matrix K;

[0088] Step S33: Core Matrix Centering

[0089] The kernel matrix is ​​centered using the following expression:

[0090]

[0091] in, For the centered approximate kernel matrix, 1 N It is an N×N matrix with all elements being 1;

[0092] Step S34: Feature decomposition and principal component extraction

[0093] Perform eigenvalue decomposition on the centered kernel matrix:

[0094]

[0095] Where V is the eigenvector matrix and Λ is the diagonal matrix;

[0096] Take the eigenvectors V corresponding to the k largest eigenvalues. k ;

[0097] The feature vector V k Projected onto principal component space:

[0098]

[0099] Where k is the number of principal components retained.

[0100] Furthermore, in step S3, the time-domain and frequency-domain features include kurtosis, peak value, and power spectral density.

[0101] Compared with the prior art, the present invention has the following advantages: This machine learning-based rolling bearing fault diagnosis method uses variational mode decomposition to reduce noise in vibration signals. To address the challenge of selecting variational mode decomposition parameters, an adaptive gray wolf algorithm is used for optimization. The permutation entropy is used to extract features from the denoised data, and a nonlinear principal component analysis algorithm is used to reduce the dimensionality of fault features. This effectively reduces the impact of noise in vibration signals on the accuracy of fault diagnosis and can significantly improve the accuracy of rolling bearing fault diagnosis. Attached Figure Description

[0102] Figure 1 This is a flowchart illustrating the rolling bearing fault diagnosis method based on machine learning in an embodiment of the present invention.

[0103] Figure 2 This is a schematic diagram of the implementation process of AGWO-VMD in an embodiment of the present invention;

[0104] Figure 3 This is a comparison chart of SNR, MSE, and RMSE before and after noise reduction using AGWO-VMD in this embodiment of the invention.

[0105] Figure 4 This is a diagram showing the diagnostic results after being diagnosed by the fault diagnosis model in an embodiment of the present invention. Detailed Implementation

[0106] The embodiments of the present invention are described in detail below. These embodiments are implemented based on the technical solution of the present invention, and provide detailed implementation methods and specific operation processes. However, the scope of protection of the present invention is not limited to the following embodiments.

[0107] Example 1

[0108] like Figure 1 As shown, this embodiment provides a technical solution: a rolling bearing fault diagnosis method based on machine learning, comprising the following steps:

[0109] Step S1: Input the vibration signals under different states into the variational mode decomposition algorithm, and use the adaptive gray wolf algorithm (AGWO) to optimize the number of modes K and the penalty factor α of the variational mode decomposition (VMD) to establish a variational mode decomposition parameter optimization model.

[0110] Step S2: Decompose the vibration signal using the optimized VMD parameters, calculate the correlation and permutation entropy of each IMF component, and filter out the effective IMF components corresponding to the original vibration signal using a double threshold.

[0111] Step S3: Input the effective IMF components into the kernelized nonlinear principal component analysis algorithm (KNPCA), take the principal components with a cumulative contribution rate of more than 90%, and construct the fault feature set of rolling bearing vibration data by combining time domain and frequency domain features.

[0112] Step S4: Train the Particle Swarm Optimized Support Vector Machine (PSO-SVM) model using the fault feature set, and use the trained model to diagnose rolling bearing faults.

[0113] In this embodiment, step S1 of the adaptive gray wolf algorithm includes the following main steps:

[0114] 1) Initial gray wolf positions are generated using Latin hypercube sampling (LHS);

[0115] 2) Calculate the initial fitness and determine the leader (α, β, δ wolf) based on the fitness;

[0116] 3) Enter the iteration loop;

[0117] 4) Output the optimal solution.

[0118] In this embodiment, the specific process of generating the initial gray wolf position using Latin hypercube sampling is as follows:

[0119] 1. Define the parameter range matrix

[0120] Construct a parameter range matrix of dimension d×2:

[0121]

[0122] Where d represents the dimension of the optimization problem (VMD requires optimizing the number of modes K and the penalty factor α, totaling 2 dimensions). Let be the minimum and maximum values ​​of the k-th dimension parameter;

[0123] 2. Generate LHS sampling matrix

[0124] 2.1 For each dimension k, generate N equally spaced intervals:

[0125]

[0126] Where i = 1, 2, ..., N are interval indices (corresponding to N individual gray wolves).

[0127] 2.2. Randomly sample one point within each interval to form a sampling sequence of dimension k:

[0128]

[0129] 2.3. Perform the above operations on all dimensions to generate d independent sampling sequences.

[0130] 3. Random permutations and combinations

[0131] 3.1 Sampling sequence S for each dimension k Perform random arrangement.

[0132] 3.2 Combine the randomly arranged sampling sequences of each dimension to generate the initial population matrix:

[0133]

[0134] Here, rand(i) represents the random permutation index of the i-th row.

[0135] In this embodiment, the formula for calculating fitness is as follows:

[0136] Fitness=ω1·MSE+ω2·(-SNR)

[0137] Where ω1 and ω2 are weighting coefficients; MSE is the mean square error between the denoised signal and the original signal; SNR is the signal-to-noise ratio of the denoised signal, which is the signal after denoising.

[0138] In this embodiment, the iterative loop includes the following steps:

[0139] 1. Dynamically adjust the convergence factor τ

[0140] Update τ using a nonlinear decay strategy:

[0141]

[0142] Where t is the iteration number; τ initial T is the initial convergence factor; max γ is the maximum number of iterations; γ is the decay coefficient.

[0143] 2. Location update

[0144] 2.1 Calculate the distance to the prey:

[0145]

[0146] 2.2 Update the location of the gray wolves:

[0147]

[0148] in, Let be the distance vector between the individual gray wolf and its prey; This represents the location of the prey, i.e., the current optimal solution. Let be the current position of the gray wolf individual at the t-th iteration, i.e., the candidate solution; It is a vector of random coefficients.

[0149] in, The calculation formula is:

[0150]

[0151]

[0152] in, [0,1] d A random vector between [variables].

[0153] 3. Adaptive mutation processing, the formula is as follows:

[0154]

[0155] Where σ is the variable asynchronous length; [0,1] d Random vectors between This represents the latest position of the prey after adaptive mutation processing, i.e., the current solution.

[0156] 4. Recalculate fitness and update α, β, and j wolves.

[0157] 5. Check the termination condition; if the condition is met, exit the loop.

[0158] It should be noted that the termination condition is reaching the maximum number of iterations or the optimal fitness being less than a threshold.

[0159] In this embodiment, in step S2, the method for filtering effective IMF components using dual thresholds is as follows:

[0160] Set a dual threshold for filtering, that is, simultaneously set a correlation coefficient threshold ρ. threshold and permutation entropy threshold PE threshold Retain those that simultaneously satisfy ρ i ≥ρ threshold And PE i ≤PE threshold The IMF component, i represents the i-th IMF component, i.e. the effective IMF component obtained by filtering.

[0161] More specifically, the formula for calculating the correlation coefficient is:

[0162]

[0163] Where x(t) is the original vibration signal with a length of N; IMF i (t) represents the i-th IMF component (modal component) obtained from VMD decomposition; It is the mean of the original signal (original vibration signal) and the i-th IMF component;

[0164] Formula for calculating permutation entropy:

[0165]

[0166] Where m is the embedding dimension; p j It represents the probability of the j-th permutation pattern.

[0167] The formula for calculating the correlation coefficient threshold is:

[0168] ρ threshold =μ ρ +0.5σ ρ

[0169] The formula for calculating the permutation entropy threshold is:

[0170] PE threshold =μ PE -1.2σ PE

[0171] Where, μ ρ μ PE σ represents the correlation coefficient and the mean of the permutation entropy. ρ σ PE For the correlation coefficient and the standard deviation of the permutation entropy, respectively.

[0172] In this embodiment, the specific steps of the kernelized nonlinear principal component analysis algorithm in step S3 are as follows:

[0173] 1. Gaussian kernel (RBF kernel) calculation:

[0174]

[0175] Where ψ is the Gaussian kernel bandwidth parameter; ||x i -x j || 2 For sample x i With x j The square of the Euclidean distance.

[0176] 2. Approximate acceleration of kernel matrix calculation:

[0177] 2.1 Random sampling: Randomly select m base samples (m << N) from N samples.

[0178] Where m is The approximate base sample size is usually taken as...

[0179] 2.2 Calculate the submatrix:

[0180] W∈R m×m : Kernel matrix between base samples, W ij =K(x) i x j ); where i, j = 1, 2, ... m;

[0181] C∈R m×(N-m) : Kernel matrix between base samples and other samples, C ij =K(x) i x j ); where i is the base sample index (there are m in total, corresponding to x1, x2…x m ); l is the index of non-base samples (a total of Nm samples, that is, the remaining samples after removing the first m base samples from the original total samples, corresponding to x). m+1 x m+2 …x N ).

[0182] 2.3 Approximate kernel matrix:

[0183]

[0184] in, This is an approximation matrix of the original kernel matrix K.

[0185] 3. Core matrix centralization:

[0186]

[0187] Among them, 1 N It is an N×N matrix consisting entirely of 1s.

[0188] 4. Characteristic Decomposition and Principal Component Extraction

[0189] Perform eigenvalue decomposition on the centered kernel matrix:

[0190]

[0191] Where V is the eigenvector matrix and Λ is the diagonal matrix.

[0192] Take the eigenvectors V corresponding to the k largest eigenvalues. k .

[0193] Project the data onto the principal component space:

[0194]

[0195] Where k is the number of principal components retained.

[0196] In this embodiment, in step S3, the time-domain and frequency-domain features include:

[0197] kurtosis K Kurtosis :

[0198] Where N is the number of sample points; x i Let i be the i-th time-domain signal sample; Standard deviation, μ is the mean.

[0199] Peak X peak :x peak =max(|x i |)

[0200] Power spectral density P:

[0201] Among them, f k For the k-th discrete frequency point; X(f k ) represents the FFT result of the signal; Δt represents the sampling interval.

[0202] Example 2

[0203] The dataset used in this embodiment is the Western Reserve University (CWRU) dataset. Acceleration data from the drive end of an SKF6205 deep groove ball bearing, sampled at a frequency of 12 kHz, with a motor speed of 1797 r / min, were selected. The bearing damage was single-point damage caused by electrical discharge machining, with fault diameters of 0.007 inches, 0.014 inches, and 0.021 inches. Gaussian white noise was used to simulate noise in the actual working environment, with a signal-to-noise ratio (SNR) of -3.

[0204] In the adaptive gray wolf algorithm, the optimization parameter dimension is set to d=2, the number of gray wolves is N=100, and the maximum number of iterations is T. max =100, initial convergence factor α initial =2. AGWO-VMD flowchart as follows Figure 2 As shown, the effect of noise reduction using AGWO-VMD is as follows: Figure 3 As shown, the signal-to-noise ratio (SNR) increased from -2.76 to 3.08, the mean square error (MSE) decreased from 0.87 to 0.23, and the root mean square error (RMSE) decreased from 0.93 to 0.48.

[0205] This embodiment uses a sliding window to construct the dataset. Each bearing fault state contains 160 samples, and each sample contains 1024 vibration acceleration data points. The first 70 sets of data are used as the training set, and the last 30 sets are used as the test set. See Table 1 for details.

[0206] Table 1. Fault Type Category Information in CWRU Dataset

[0207] Fault type Fault diameter (inch) Type markup normal — 1 Inner ring fault 0.007 2 Inner ring fault 0.014 3 Inner ring fault 0.021 4 Outer ring fault 0.007 5 Outer ring fault 0.014 6 Outer ring fault 0.021 7 Rolling element failure 0.007 8 Rolling element failure 0.014 9 Rolling element failure 0.021 10

[0208] Diagnostic results as follows Figure 4 As shown, the horizontal axis represents the original label of the fault type, and the vertical axis represents the label classified after model diagnosis. It can be seen that all labels are correctly classified, with an accuracy rate of 100%.

[0209] In summary, the machine learning-based rolling bearing fault diagnosis method described in the above embodiments reduces noise in vibration signals by utilizing variational mode decomposition. To address the challenge of selecting variational mode decomposition parameters, an adaptive gray wolf algorithm is employed for optimization. Permutation entropy is used to extract features from the denoised data, and nonlinear principal component analysis is applied to reduce the dimensionality of fault features. These methods effectively reduce the impact of noise in vibration signals on fault diagnosis accuracy and significantly improve the accuracy of rolling bearing fault diagnosis.

[0210] Although embodiments of the present invention have been shown and described above, it is understood that the above embodiments are exemplary and should not be construed as limiting the present invention. Those skilled in the art can make changes, modifications, substitutions and variations to the above embodiments within the scope of the present invention.

Claims

1. A machine learning-based method for diagnosing rolling bearing faults, characterized in that, Includes the following steps: Step S1: VMD parameter optimization Vibration signals under different states are input into the variational mode decomposition algorithm, and the adaptive gray wolf algorithm is used to optimize the variational mode decomposition parameters, which include the number of modes K and the penalty factor α. Step S2: Screening of Effective IMF Components The vibration signal is decomposed using the optimized variational mode decomposition parameters to obtain multiple IMF components. The correlation coefficient and permutation entropy of each IMF component are calculated, and the effective IMF components corresponding to the original vibration signal are selected by double thresholding. Step S3: Feature Extraction and Dimensionality Reduction The effective IMF components are input into the kernelized nonlinear principal component analysis algorithm. The principal components whose cumulative contribution rate exceeds the set threshold are selected, and the fault feature set of rolling bearing vibration data is constructed by combining time domain and frequency domain features. Step S4: Fault Diagnosis The particle swarm optimization support vector machine model is trained using the fault feature set, and the trained model is used for fault diagnosis of rolling bearings. In step S1, the specific processing procedure of the adaptive gray wolf algorithm is as follows: Step S11: Generate the initial gray wolf position using Latin hypercube sampling; Step S12: Calculate the initial fitness and determine the leader based on the fitness, where the leader is α, β, or δ wolf; Step S13: Enter the iteration loop; Step S14: Output the optimal solution; In step S3, the specific processing procedure of the kernelized nonlinear principal component analysis algorithm is as follows: Step S31: Gaussian kernel calculation The Gaussian kernel K(x) is calculated using the following formula. i x j ): Where ψ is the Gaussian kernel bandwidth parameter; ||x i -x j || 2 For sample x i With x j The square of the Euclidean distance; Step S32: Approximate accelerated kernel matrix calculation Step S321: Random sampling Randomly select m base samples from N samples, where m << N and m is the approximate number of base samples; Step S322: Calculate the submatrix Calculate the kernel matrix between the base samples: W ij =K(x i ,x j ) Where W∈R m×m i, j = 1, 2, ... m; Calculate the kernel matrix between the base sample and other samples: C il =K(x i ,x m+l ) Where C∈R m×(N-m) , where i is the base sample index; l is the non-base sample index; Step S323: Approximate kernel matrix The kernel matrix is ​​approximated using the following expression: in, This is an approximation matrix of the original kernel matrix K; Step S33: Core Matrix Centering The kernel matrix is ​​centered using the following expression: in, For the centered approximate kernel matrix, 1 N It is an N×N matrix with all elements being 1; Step S34: Feature decomposition and principal component extraction Perform eigenvalue decomposition on the centered kernel matrix: Where V is the eigenvector matrix and Λ is the diagonal matrix; Take the eigenvectors V corresponding to the k largest eigenvalues. k ; The feature vector V k Projected onto principal component space: Where k is the number of principal components retained.

2. The rolling bearing fault diagnosis method based on machine learning according to claim 1, characterized in that, In step S11, the specific processing procedure is as follows: Step S111: Define the parameter range matrix Construct a parameter range matrix of dimension d×2: Where d represents the dimension of the optimization problem. Let K be the minimum and maximum values ​​of the k-th dimension parameter. There are two parameters to be optimized: the number of modes K and the penalty factor α. The value of d is 2. Step S112: Generate LHS sampling matrix For each dimension k, generate N equally spaced intervals: Where i = 1, 2, ..., N are interval indices, corresponding to N individual gray wolves; Randomly sample a point within each interval to form a sampling sequence of dimension k: Perform the above operation on all dimensions to generate d independent sampling sequences; Step S113: Random permutation and combination Sampling sequence S for each dimension k Randomly arrange the sampled sequences of each dimension and combine them to generate the initial population matrix: Here, rand(i) represents the random permutation index of the i-th row.

3. The rolling bearing fault diagnosis method based on machine learning according to claim 1, characterized in that, In step S12, the fitness calculation formula is as follows: Fitness=ω1·MSE+ω2·(-SNR) Where ω1 and ω2 are weighting coefficients; MSE is the mean square error between the denoised signal and the original vibration signal; and SNR is the signal-to-noise ratio of the denoised signal.

4. The rolling bearing fault diagnosis method based on machine learning according to claim 1, characterized in that, In step S13, the specific processing procedure is as follows: Step S131: Dynamically adjust the convergence factor τ Update τ using a nonlinear decay strategy: Where t is the iteration number; τ initial T is the initial convergence factor; max γ is the maximum number of iterations; γ is the decay coefficient. Step S132: Location Update Calculate the distance to the prey: Updated Grey Wolf location: in, Let be the distance vector between the individual gray wolf and its prey; This represents the location of the prey, i.e., the current optimal solution. Let be the current position of the gray wolf individual at the t-th iteration, i.e., the candidate solution; The coefficient vector is calculated using the following formula: in, [0,1] d Random vectors between; Step S133: Perform adaptive mutation processing, as shown in the following formula: Where σ is the variable asynchronous length; [0,1] d Random vectors between This represents the latest position of the prey after adaptive mutation processing, i.e., the current solution. Step S134: Recalculate fitness and update α, β, and δ wolves; Step S135: Check the termination condition. If the termination condition is met, exit the loop. The termination condition is reaching the maximum number of iterations or the optimal fitness is less than the threshold.

5. The rolling bearing fault diagnosis method based on machine learning according to claim 1, characterized in that, In step S2, a dual threshold screening is set, that is, the correlation coefficient threshold ρ is set simultaneously. threshold and permutation entropy threshold PE threshold Retain those that simultaneously satisfy ρ i ≥ρ threshold And PE i ≤PE threshold The IMF component, i represents the i-th IMF component, i.e. the effective IMF component obtained by filtering.

6. The rolling bearing fault diagnosis method based on machine learning according to claim 5, characterized in that, The formula for calculating the correlation coefficient threshold is: r threshold =μ ρ +0.5σ ρ The formula for calculating the permutation entropy threshold is: PE threshold =μ PE -1.2σ PE Where, μ ρ μ PE σ represents the correlation coefficient and the mean of the permutation entropy. ρ σ PE For the correlation coefficient and the standard deviation of the permutation entropy, respectively.

7. The rolling bearing fault diagnosis method based on machine learning according to claim 6, characterized in that, In step S2, the correlation coefficient is calculated using the following formula: Where x(t) is the original vibration signal with a length of N; IMF i (t) represents the i-th IMF component obtained from VMD decomposition; The mean of the original vibration signal and the i-th IMF component; The formula for calculating permutation entropy is: Where m is the embedding dimension; p j Let be the probability of the j-th permutation pattern.

8. The rolling bearing fault diagnosis method based on machine learning according to claim 1, characterized in that, In step S3, the time-domain and frequency-domain features include kurtosis, peak value, and power spectral density.

Images