Wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration

Abstract

The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration comprises the following steps: step one: analyzing the characteristics of wind turbine SCADA data, preprocessing the wind turbine SCADA data, constructing a multivariate time series sample set containing secondary classification labels based on the matrix profile analysis method of time series; step two: a stacked LSTM model is used to design a fault diagnosis model suitable for the multivariate time series sample set, which can diagnose the overall multivariate time series sample and preliminarily predict the state of each time sample; step three: a hierarchical fault probability calibration strategy is used to fuse multivariate information to correct the probability of each time sample belonging to different fault types, and the class corresponding to the maximum probability value after correction is selected as the final fault diagnosis result. The method can fully fuse multivariate information to improve the accuracy and reliability of fault diagnosis.

Description

Technical Field

[0001] This invention belongs to the field of wind turbine fault diagnosis, specifically involving a wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration. Background Technology

[0002] Most current wind turbines are installed in remote areas with poor transportation or offshore. Due to the harsh operating environment and complex operating procedures, the failure rate of wind turbines is often higher than that of traditional thermal and hydroelectric turbines. This affects the power generation efficiency of wind turbines and increases the operation and maintenance costs of wind farms. Therefore, it is necessary to adopt reliable methods to accurately and promptly identify faults in wind turbines, thereby reducing downtime for maintenance.

[0003] Most current wind farms are equipped with Supervisory Control and Data Acquisition (SCADA) systems to monitor the operating status of wind turbines in real time, providing abundant data resources for data-driven wind turbine fault diagnosis. Compared to wind turbine fault diagnosis methods based on physical models and signal processing, data-driven methods use machine learning and deep learning to diagnose faults, offering advantages such as broad applicability, high accuracy, and no need for additional investment, thus showing broad application prospects. However, most current data-driven wind turbine fault diagnosis methods identify fault types by mining the mapping relationships between variables from the perspective of data correlation, and most methods ignore the inherent time-series attributes of wind turbine condition monitoring quantities. Even some methods that consider the inherent short-term dependencies of variables from a time-series perspective, such as recurrent neural networks and long short-term memory networks (LSTM), only consider the short-term dependencies of single variables on the time axis, ignoring the correlations between multiple time-series variables at the same time point. This results in some effective information not being fully extracted.

[0004] In addition, most current mainstream fault diagnosis methods based on classification algorithms determine the fault type by the probability value output by the model. However, due to measurement errors during data acquisition, random variations in wind speed, and other uncertain environmental factors, the quality of data collected by SCADA systems is often difficult to guarantee, which seriously affects the reliability of the model's fault diagnosis.

[0005] Furthermore, wind turbine failures are low-probability, sporadic events, resulting in a limited number of relevant failure samples with severe imbalances. This makes it difficult for trained fault diagnosis models to meet performance requirements. While some scholars have proposed resampling data to alleviate the imbalance, these methods are mostly designed for single-fault diagnosis, and their accuracy in diagnosing multiple types of faults in wind turbines still needs improvement. Summary of the Invention

[0006] To address the challenges in wind turbine fault diagnosis, this invention proposes a wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration, using SCADA data of the wind turbine as a foundation. This method can fully integrate multivariate source information to improve the accuracy and reliability of fault diagnosis.

[0007] The technical solution adopted in this invention is as follows:

[0008] A wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration includes the following steps:

[0009] Step 1: Analyze the characteristics of wind turbine SCADA data, preprocess the wind turbine SCADA data, and construct a multivariate time series sample set with secondary classification labels based on the matrix profile analysis method of time series.

[0010] Step 2: A fault diagnosis model suitable for multivariate time series sample sets was designed using a stacked LSTM model. This model performs fault diagnosis on the overall multivariate time series sample set and makes a preliminary diagnosis on the sample status at each time step within the set.

[0011] Step 3: Employ a hierarchical fault probability correction strategy to integrate multivariate information to correct the probability of samples belonging to different fault types at each time point, and select the class corresponding to the highest corrected probability value as the final fault diagnosis result.

[0012] Step one includes the following steps:

[0013] Step 1.1: Preprocess the SCADA data of the wind turbine to eliminate outliers and missing values;

[0014] Step 1.2: Combining the SCADA data and fault logs of the wind turbine, label the SCADA sampled data at each preprocessed moment; where no time is recorded in the fault log, it indicates normal operation and is labeled with the integer 0; other faults recorded in the fault log are labeled with different integers, and the corresponding sample label set is Y = [0, 1, 2, ..., K], where K is the number of fault categories.

[0015] Step 1.3: Since the SCADA system of wind turbine units typically collects hundreds of monitoring data, but not all of them are meaningful for fault diagnosis, it is necessary to combine domain knowledge with correlation analysis to select strongly correlated monitoring data as feature variables. This will reduce the amount of data to be analyzed while retaining as much effective information as possible. The selected strongly correlated monitoring data will then serve as the basic variables for constructing the sample set for multivariate time series analysis, as detailed below:

[0016] Pearson correlation analysis is used to select characteristic variables, and the formula is as follows:

[0017]

[0018] In the formula: A and B represent the sets of candidate variables and the target variable, respectively; a i b i Let μ represent the i-th element in A and B respectively; A and μ B P represents the mean of the two feature subsets, respectively. r P represents the Pearson correlation coefficient between two feature subsets; r The value of reflects the degree of correlation between two feature subsets. Its value ranges from [-1, 1]. A positive value indicates a positive correlation, a negative value indicates a negative correlation, and a value of 0 indicates no correlation. The closer the absolute value is to 1, the stronger the correlation between the two feature subsets.

[0019] Step 1.4: Standardize the selected monitoring variables, including wind speed, wind direction, various rotational speeds, various temperature monitoring quantities, electrical monitoring quantities, and output power, using the maximum and minimum value scaling method shown in Equation (2);

[0020]

[0021] In equation (2): x is the data to be processed, x * For the processed data, x max and x min These are the maximum and minimum values ​​of the monitored data for this variable, respectively.

[0022] Step 1.5: Based on the selected feature variable data, a sliding window with a fixed window width is used to sample at a certain step size to obtain a multivariate time series data matrix X containing multiple consecutive time points. i Furthermore, with X i The corresponding labels at each time point are used to construct the label vector Y. i As X i The first-level label, and through Y i Perform an OR operation on the labels inside to get X iThe second-level labels, then the newly constructed multivariate time series sample set takes the following form:

[0023]

[0024]

[0025] In equation (3): d represents the number of selected feature variables; L represents the length of the sampling window; N represents the total number of multivariate time series samples;

[0026] Vector Y i The element in the f array takes a value from [0, 1, 2, ..., K], where 0, 1, 2, ..., m represent the labels set according to the fault log, 0 represents no fault, and K represents the label for the Kth type of fault. i represents the secondary label of the multivariate fault sample, with a value of 0 indicating that there is no fault within the sliding window sampling period, and 1 indicating that there is a fault within the sliding window sampling period. i represents the i-th multivariate time series matrix sample constructed.

[0027] To facilitate differentiation, the original SCADA system samples data at each moment as a SCADA sample. Each SCADA sample contains the values ​​of multiple monitoring quantities collected by different sensors at the same moment, and corresponds to a primary label. Each multivariate time series sample contains L SCADA samples and their labels, as well as corresponding secondary labels.

[0028] Step 1.6: Divide the constructed multivariate time series into training set, validation set and test set according to a certain ratio. The training set is used for model training and parameter tuning, the validation set is used for performance evaluation of the trained model, and the test set contains several samples with known fault types for verifying the actual effect of the model.

[0029] In step 1.1, the 3δ criterion is used to initially identify outliers caused by sampling anomalies in various monitoring quantities, and interpolation is used to replace and supplement outliers and missing values ​​in the original data.

[0030] The 3δ criterion, also known as the Raida criterion, is a statistical method used to detect outliers. This criterion assumes that the data follows a normal distribution. Since the probability of standard normal distribution data falling within the interval (μ-3σ, μ+3σ) is 99.74%, the probability of data falling outside this range is less than 0.3%. Therefore, the 3δ criterion considers data points outside the (μ-3σ, μ+3σ) interval to be outliers, requiring investigation and handling. Here, μ represents the mean of the data set, and σ is the standard deviation, reflecting the dispersion of the data.

[0031] Interpolation is a commonly used numerical analysis method used to estimate function values ​​at unknown points by constructing an interpolation function when data points are known. Common interpolation methods include Lagrange interpolation and Newton interpolation. This invention uses the commonly used Lagrange interpolation method to handle outliers and missing values. Its principle is as follows:

[0032] Suppose we have n data points (x1, y1), (x2, y2), ..., (x... n ,y n The task is to construct an (n-1)th degree polynomial P(x) on these data points, such that P(x) i )=y i (i = 1, 2, ..., n). Then the Lagrange interpolation formula for P(x) is:

[0033]

[0034] Among them, y i x represents the function value of the i-th data point; i Let represent the value of the independent variable at the i-th data point; the π symbol in the formula indicates multiplication over all j ≠ i. Handling outliers and missing values ​​can improve data quality and thus enhance the performance of subsequent model training.

[0035] Step two includes the following steps:

[0036] Step 2.1: Based on the stacked LSTM model, mine the dependencies between time series data and faults across multiple time periods. The input to the stacked LSTM model is the data matrix X from the multivariate time series samples. i As input, the first-level label vector Y i and second-level tag f i These are the output training models. The loss function corresponding to the first-level labels is the multivariate cross-entropy loss function L. m Second-level tag f i The corresponding loss function is the binary cross-entropy loss function L. b The corresponding expressions are shown in equations (5) and (6) respectively:

[0037]

[0038]

[0039] In equations (5) and (6): N vThe symbols represent the number of multivariate time series samples; L represents the number of SCADA samples contained in each multivariate time series sample; the superscript c indicates that the corresponding classification label is c; K represents the total number of labels contained in the first-level label; the subscript i represents the i-th multivariate time series sample; and the subscript t represents the sample collected by the SCADA system at time t in each multivariate time series sample, meaning that each sample corresponds to a unique t. This notation is used because when constructing multivariate time series samples using sliding window sampling, SCADA samples at the same time may be assigned to different multivariate time series samples, but the unique timestamp t can still accurately identify the sample.

[0040] This indicates whether the label of the t-th SCADA sample in the i-th multivariate time series sample is c. This represents the probability that the model classifies the label of the t-th SCADA sample in the i-th multivariate time series sample as c. Log() is used to calculate the logarithm. i and p i Let represent the label of the i-th sample as positive and its corresponding probability when classified by secondary labels;

[0041] Stacked LSTM network is a deep learning model used for processing time series data. It enhances the model's ability to mine data by stacking multiple LSTM layers. The schematic diagram of the stacked network is shown below. Figure 3 As shown, the stacked LSTM model for top-level fault diagnosis is used to make a preliminary prediction of the overall fault probability of the multivariate time series samples; while the stacked LSTM model for internal fault probability diagnosis is used to predict the fault probability of samples at each time step within the multivariate time series samples. The number of stacked layers in these two stacked LSTM models are N1 and N2, respectively, and their values ​​need to be determined based on the actual performance of the two models.

[0042] Step 2.2: Optimize the parameters of the stacked LSTM model using a grid search method and appropriate evaluation metrics. Evaluation metrics include the accuracy metric A based on multivariate time series samples. fir The average diagnostic accuracy index A for all time series samples at each time point. sec Average false positive rate index F NR Average false alarm rate index F PR ;

[0043]

[0044]

[0045]

[0046]

[0047] In the formula: L represents the length of the sliding window when constructing the multivariate time series samples, which includes L original samples collected by the SCADA system; N TP N represents the number of normal samples identified as normal. TN This represents the number of fault classes identified as fault classes in a multivariate time series sample. This represents the number of fault classes that are incorrectly identified as normal classes in all time-series samples within the i-th multivariate time series sample; This represents the number of normal samples that are incorrectly identified as faulty samples in all time steps within the i-th multivariate time series sample.

[0048] Of the above indicators: A fir This indicates the classification accuracy of the stacked LSTM model for a multivariate time series sample population; A sec The average false positive rate (F) represents the accuracy of the model in classifying samples across all time points within a multivariate time series sample. NR Used to reflect the missed detection rate of the model in fault diagnosis at each time step; average false alarm rate index F PR This is used to reflect the false alarms in the model's fault diagnosis at each time point.

[0049] Step three includes the following steps:

[0050] Step 3.1: Calculate the global failure probability P of the multivariate time samples. i As top-level prior information, the probability of individual failures The values ​​are corrected to eliminate the influence of individual extreme cases on the diagnostic results; an empirical weighting method is used to... The value is corrected once:

[0051]

[0052] In equation (11): ω represents the probability of belonging to the c-th type of fault at time t within the i-th multivariate time series sample after the first correction; ω is the probability fusion weight, and the value of ω reflects the confidence of the model in the individual prediction results; To determine the expected values ​​of the probabilities of various types of faults under the condition of global sample faults in a multivariate time series; Indicates taking The superscript c with the highest probability value means that the maximum probability of fault classification for sample at time t is consistent with the global sample fault probability value.

[0053] The fault classification probability of the sample at time t after one correction using the top-level probability is: Let represent the probabilities of the sample sampled at time t within the i-th multivariate time series sample after one correction, belonging to the fault (0, 1, 2, ..., K). The corrected probability values ​​take into account the global information statistical distribution of the multivariate time series samples.

[0054] Step 3.2: Fault probability estimation method based on statistical inference, including the following steps:

[0055] Step 3.2.1: Based directly on the preprocessed training data, i.e., the SCADA sampling data arranged in the original time series order, construct K different training sets D according to the categories. c c = 0, 1, 2, ..., K. The first d columns of the dataset represent different monitoring variables, and the last column is the corresponding fault sample label. Since the number of samples of each type is not necessarily the same, D0 to D... K The number of rows in the matrices is different, as shown in equation (12): D0 and D c The subscripts l and m respectively indicate that there are l normal type samples and m c type fault samples;

[0056]

[0057] In equation (12): This represents the standardized value of the d-th monitoring variable in the m-th sample of the c-th type of fault. Step 3.2.2: A non-parametric kernel density estimation algorithm is used to learn from historical data of different fault types, to mine the conditional dependencies between various monitoring data and different fault labels, and to obtain the conditional probability density function f for each type of fault. c (x), then given a new sample conditions, Let f represent d different monitored variables collected by the SCADA system at time t. c (x) can estimate the probability that the sample belongs to different types of faults, and the probability is expressed as:

[0058]

[0059] In equation (13): This represents the monitored quantity x at a given input time t. t In this case, the conditional probability density function f of type c faults is used. c (x) estimates the probability that it belongs to type c fault; m is the number of training samples; h is a hyperparameter, called bandwidth or window; d represents the number of monitored variables; φ(x) is the Gaussian kernel function; D c The training set is constructed from all sample data of the c-th type of fault; h d The meaning is the multiplication of d h's; x tLet x represent the input condition vector at a given time t. i D represents c The condition vector acquired at the i-th time step.

[0060] By estimating the probability density functions for different fault types under the same input conditions, we can obtain the probabilities belonging to different fault types. Let q represent the probabilities that a sample at time t belongs to the 0th, 1st, 2nd, ..., Kth time step, respectively, based on kernel density estimation. t It estimates the probability using all samples of each type of fault in the training sample space. The estimation result contains the true distribution information of each type of fault sample and the individualized information of each sample. Therefore, this fault probability can be used as the underlying prior probability to further refine p. t,i Make corrections.

[0061] Step 3.3: Fault probability calibration based on Bayes' theorem, in Based on the conditional kernel density estimation results of various faults, q t As the underlying prior probability, the probability of a sample belonging to different faults at time t is further corrected. The specific method is as follows:

[0062]

[0063] In equation (14): Let be the posterior probability of the i-th multivariate time series sample belonging to the c-th type of fault at time t, after secondary correction. This represents the underlying prior probability that a sample at time t belongs to fault c, obtained based on different kernel density estimates; and the posterior probability set of a sample at time t belonging to each class, obtained from these values. back, Let and represent the probabilities that the sample at time t belongs to the 0th, 1st, 2nd, ..., Kth time groups after secondary correction. Finally, the sample with the highest probability value is selected as the diagnosed fault type, and the corresponding reliability probability is .

[0064] This invention provides a wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration, with the following technical advantages:

[0065] 1) The matrix contour analysis method based on time series provided by this invention constructs multivariate time series samples with secondary labels, which can effectively alleviate the impact of noise and missing values ​​in the original data on model training.

[0066] 2) This invention utilizes a stacked LSTM model trained with multivariate time series samples, which can not only effectively learn the short-term dependencies within the data, but also extract the correlation between different variables at the same time point, effectively improving the depth of information mining within the data.

[0067] 3) This invention proposes a hierarchical fault probability correction strategy that combines Bayesian theory, which can fully integrate multi-source information to improve the accuracy and reliability of fault diagnosis. Attached Figure Description

[0068] Figure 1 A flowchart for fault diagnosis of wind turbine units combining LSTM and hierarchical Bayesian inference.

[0069] Figure 2 A schematic diagram for constructing multivariate time series samples.

[0070] Figure 3 This is a general framework diagram for fault diagnosis and probabilistic correction based on stacked LSTM.

[0071] Figure 4 Classification results of multivariate time series samples based on secondary labels

[0072] Figure 5 This is a diagnostic result of the internal sample status of the same multivariate time series before (left) and after (right) probability correction. Detailed Implementation

[0073] This invention proposes a fault diagnosis method based on stacked LSTM and hierarchical probability calibration, using SCADA data from wind turbine generators. The method comprises three steps to diagnose multiple types of faults in wind turbine generators. First, the characteristics of the SCADA data are analyzed, and the data is preprocessed using professional knowledge and correlation analysis. Based on this, matrix profile analysis and sliding window sampling are combined to construct a multivariate time-series sample set with secondary classification labels. Second, a fault diagnosis model suitable for the multivariate time-series sample set is designed using a stacked LSTM model, which can perform fault diagnosis on the overall multivariate time-series sample set while simultaneously making preliminary predictions of the sample status at each time step. Finally, hierarchical probability calibration is used to fuse multivariate information to correct the probability of samples at each time step belonging to different fault types, and the class corresponding to the highest corrected probability value is selected as the final fault diagnosis result. The basic flowchart of the proposed method is shown below. Figure 1 As shown.

[0074] This invention can be divided into three steps: data preprocessing and construction of a multivariate time series sample set, fault diagnosis based on stacked LSTM, and fault probability hierarchical calibration based on empirical weighting and Bayes' theorem. The detailed technical steps are as follows:

[0075] (I): Data preprocessing and construction of multivariate time series sample sets:

[0076] 1) Due to various factors, the quality of raw data collected by wind turbine SCADA systems is usually poor, requiring preprocessing to eliminate outliers and missing values. This invention employs the 3δ criterion to initially identify outliers caused by sampling anomalies, and uses interpolation to replace and supplement outliers and missing values ​​in the raw data.

[0077] 2): Standardize each type of monitoring variable by using maximum and minimum value scaling as shown in equation (2), where x is the data to be processed, x * For the processed data, x max and x min These are the maximum and minimum values ​​of the monitored data for this variable, respectively.

[0078]

[0079] 3) Combining SCADA data and wind turbine fault logs, the fault logs record the fault alarm information of the wind turbines. Labels are set for the SCADA sampled data at each preprocessed moment. Moments not recorded in the fault logs indicate normal operation and are labeled with the integer 0. Other faults recorded in the fault logs are labeled with different integers. The corresponding sample label set is Y = [0, 1, 2, ..., K], where K is the number of fault categories.

[0080] 4): In order to improve the learning efficiency of the model and avoid the influence of irrelevant variables on the model learning, the mapping relationship between variables and faults is mined from the sample data. Domain knowledge and correlation analysis are combined to select feature variables to construct a sample set.

[0081] 5): Based on the selected feature variable data, a sliding window with a fixed window width is used to sample at a certain step size to obtain a multivariate time series data matrix X containing multiple consecutive time points. i Furthermore, with X i The corresponding labels at each time point are used to construct the label vector Y. i As X i The first-level label, and through Y i Perform an OR operation on the labels inside to get X i The second-level labels, then the newly constructed multivariate time series sample set takes the following form:

[0082]

[0083] In equation (3): d represents the number of selected feature variables, L represents the length of the sampling window, N represents the total number of multivariate time series samples, and vector Y... i The range of values ​​for the elements in f is [0, 1, 2, ..., m]. i This serves as the secondary label for the multivariate fault samples, with a value of 0 indicating no faults within the sliding window sampling period and 1 indicating the presence of faults within the sliding window sampling period. The specific multivariate time series sampling construction process is as follows... Figure 2 As shown, the sliding window moves a distance of v each time, which can be set according to requirements.

[0084] 6): Divide the constructed multivariate time series into training set, validation set and test set according to a certain ratio. The training set is used for model training and parameter tuning, the validation set is used for performance evaluation of the trained model, and the test set contains several samples with known fault types for verifying the actual effect of the model.

[0085] (II): Fault probability diagnosis based on stacked LSTM:

[0086] 1) Based on a stacked LSTM model, mine the correlation between time series data and faults across multiple time periods. The model input is a data matrix X from a multivariate time series sample. i As input, the first-level label vector Y i and second-level tag f i These are the output training models. The loss function corresponding to the first-level labels is the multivariate cross-entropy loss function L. m The loss function corresponding to the second-level labels is the binary cross-entropy loss function L. b The corresponding expressions are shown in equations (5) and (6):

[0087]

[0088]

[0089] In equations (5) and (6): N v The meanings of L and K are: the number of multivariate time series samples, L represents the number of SCADA samples in each multivariate time series sample, K represents the total number of labels in the first-level label, subscript i represents the i-th multivariate event sequence sample, subscript t represents the t-th SCADA sample in each multivariate time series sample, and superscript c represents the corresponding label c. This indicates whether the label of the t-th SCADA sample in the i-th multivariate time series sample is c. This represents the probability that the model predicts the label of the t-th SCADA sample in the i-th multivariate time series sample to be c. Log() is the logarithmic calculation. It should be noted that the two parts of the output of the stacked LSTM model contain different physical meanings. The output based on the secondary label reflects whether there is a fault alarm within the multivariate time series sample, while the output based on the primary label reflects the probability of the fault alarm belonging to a certain category.

[0090] 2) Furthermore, based on experience, a set of values ​​for various parameters of the stacked LSTM model, including the number of hidden units, the number of LSTM layers, and other hyperparameters (such as learning rate and batch size), is set. Using the validation set as a basis, a grid search method and corresponding evaluation metrics are employed to optimize the parameters of the stacked LSTM model. The model evaluation metrics adopted by the designed method include the accuracy metric A based on multivariate time series samples. fir The average diagnostic accuracy index A for all time series samples at each time point. sec Average false positive rate index F NR Average false alarm rate index F PR :

[0091]

[0092]

[0093]

[0094]

[0095] Where: N TP N represents the number of normal samples identified as normal. TN This represents the number of fault classes identified as fault classes in a multivariate time series sample. This represents the number of fault classes that are incorrectly identified as normal classes in all time-series samples within the i-th multivariate time series sample. This represents the number of normal samples incorrectly identified as faulty samples across all time points within the i-th multivariate time series sample. In the above metric, A... fir A represents the classification accuracy of the stacked LSTM model for multivariate time series samples. sec The average false positive rate (F) represents the accuracy of the model in classifying samples at each time step. NR The average false alarm rate (F) is used to reflect the missed detection rate of the model in fault diagnosis at each time step. PR This is used to reflect the false alarms in the model's fault diagnosis at each time point.

[0096] The above evaluation of the LSTM model focuses more on assessing the accuracy of secondary label diagnosis. This is mainly because assessing the accuracy of primary label diagnosis is a simple binary classification problem with a large amount of sample data, theoretically leading to better training and classification accuracy. Secondary label classification, however, is a multi-class problem. The quality of the training data and the imbalance in the number of fault samples across different categories significantly impact model performance and therefore deserve greater attention.

[0097] (III): Fault Probability Hierarchical Calibration Based on Empirical Weighting and Bayes' Theorem

[0098] As mentioned earlier, the stacked LSTM model outputs the fault probability P of multivariate time series samples through different loss functions and fully connected layers. i and the fault classification probability at each time point within this multivariate time series. t represents the sampling time, t = 1, 2, 3, ..., N all N all To verify the total number of records included in the set, N all =N v *L. However, due to factors such as training data quality and model parameter settings, the prediction results of deep learning models may not accurately reflect the true distribution of samples, i.e., overfitting and underfitting may occur. Therefore, it is necessary to combine multivariate information to correct the model's predicted failure probability. This invention designs a two-stage failure probability correction method as follows:

[0099] ①: Due to the multivariate time sample probability P estimated by stacked LSTM i Relative to the failure probability p at each internal time point t,i It belongs to the top-level probability, theoretically p t,i Internal elements It should follow the fault distribution as a global sample, therefore the global fault probability P of the multivariate time sample can be... i As top-level prior information, the probability of individual failures The values ​​are corrected to eliminate the influence of extreme individual cases on the diagnostic results. Here, an empirical weighting method is used to correct the bias. The value is corrected once:

[0100]

[0101] In the formula: The probability of a fault belonging to type c at time t within the i-th multivariate time series sample is the result after the first correction. To determine the expected values ​​of the probabilities of various types of faults under the condition of global sample faults in a multivariate time series. Indicates taking The superscript 'c' represents the highest probability value, ensuring that the maximum fault classification probability of the sample at time t is consistent with the global sample fault probability value. ω is the probability fusion weight; the value of ω reflects the model's confidence in the individual prediction results. Here, ω is set to 0.8 based on experience. The fault classification probability of the sample at time t after one correction using the top-level probability is... The corrected probability value considers both the global statistical distribution of information from multivariate time samples and, compared to p... t,i It has better reliability.

[0102] ②: To more reliably obtain the true distribution of various faults and eliminate the risks of overfitting and underfitting introduced during deep learning, this invention also estimates the fault probability of the samples from a statistical inference perspective, using this as the underlying prior probability to perform a secondary correction on the fault classification probability. The fault probability estimation method based on statistical inference is as follows:

[0103] Instead of using matrix profiling techniques to construct multivariate time series samples, this method directly constructs K different training sets D based on the original training data, categorized by class. c (c = 0, 1, 2, ..., K), the first d columns of the dataset represent different monitoring variables, and the last column is the corresponding fault sample label. Since the number of samples for each type is not necessarily the same, D0~D K The number of rows in the matrices are different, as shown in equation (12) D0 and D c The subscripts l and m indicate that there are l normal type samples and m c type fault samples, respectively.

[0104]

[0105] In the formula: This represents the standardized value of the d-th monitoring variable in the m-th sample of the c-th type of fault.

[0106] Furthermore, a nonparametric kernel density estimation algorithm is employed to learn from historical data of different fault types, mining the conditional dependencies between various monitoring data and different fault labels to obtain the conditional probability density function f for each type of fault. c (x), then given a new sample Under different conditions, using different f c (x) can estimate the probability that the sample belongs to different types of faults, and the probability is expressed as:

[0107]

[0108] In the formula: This represents the monitored quantity x at a given input time t. tIn this case, the probability of a fault belonging to class c is estimated using the conditional probability density function fc(x) of class c faults; m is the number of training samples; h is a hyperparameter, called bandwidth or window; d has the same meaning as in the data matrix above, representing the number of monitored variables; φ(x) is the kernel function, and this invention uses the Gaussian kernel function. By estimating the probability density functions of different fault types under the same input conditions, the probabilities belonging to different fault types can be obtained. q t It estimates the probability using all samples of each type of fault in the training sample space. The estimation result contains the true distribution information of each type of fault sample and the individualized information of each sample. Therefore, this fault probability can be used as the underlying prior probability to further refine p. t,i Make corrections. It should be noted that q... t The global distribution information of the overall sample space was not considered (because the distribution of various types of samples in the overall sample space is unbalanced), therefore, it was used alone. Fault diagnosis carries a high risk of overfitting, and the accuracy of simple kernel density estimation methods is also affected by the number of fault samples of various types. Therefore, it should be used with caution. Perform fault diagnosis.

[0109] ③: Utilizing the top-level prior probability P i Obtained by performing a correction While considering the influence of the statistical distribution of the sample space can effectively eliminate the interference of extreme values ​​on fault diagnosis, it cannot effectively describe the individualized information of the samples. Therefore, this invention further utilizes Bayes' theorem... Based on the conditional kernel density estimation results of various faults, q t The probability of a sample belonging to different faults at time t is corrected a second time using the underlying prior probability. The specific method is as follows:

[0110]

[0111] In the formula: Let be the posterior probability, after secondary correction, that the i-th multivariate time series sample belongs to the c-th type of fault at time t. This represents the underlying prior probability that a sample at time t belongs to fault c, obtained based on different kernel density estimates. Once the posterior probability set of a sample at time t belonging to each class is obtained... Then, the fault type with the highest probability value is selected as the fault type to be diagnosed, and the corresponding reliability probability is...

[0112] Verification Example:

[0113] (1) Validation data description:

[0114] To verify the effectiveness of the proposed method, data collected from the SCADA system of a 1.5MW wind turbine at an onshore wind farm from November 1, 2017 to October 31, 2018 was used for reliability analysis. The data sampling interval was 5 minutes, and information on 113 variables, including environmental variables and turbine status variables, was monitored, resulting in 101,120 valid data records. Although the wind turbine did not experience any prolonged downtime for maintenance due to serious faults during the year, the SCADA system's fault logs and manual maintenance records showed that the turbine still frequently experienced minor but frequent faults such as temperature exceeding limits, voltage or current exceeding limits, or drive abnormalities. According to the manufacturer's fault handling manual for this type of wind turbine, the manufacturer classifies these faults into 5 major categories and 29 subcategories. For ease of study, the major categories are used as fault classification labels in the example. The number of fault alarms for each major category and the subcategories of faults included in each major category during the year are shown in Table 1.

[0115] Table 1. Classification of Main Fault Types and Statistics of Fault Frequency for Wind Turbine Units

[0116]

[0117] Using the multivariate time series sample construction method described earlier, the sliding window length L = 144 was chosen in the example, meaning each multivariate time series window contains 12 consecutive hours of data, with each sliding distance being 12 sampling points, or 1 hour. Combining professional knowledge and correlation analysis, 25 variables, including wind speed, power, main shaft speed, gearbox oil temperature, generator winding temperature, hydraulic pressure, and pitch position, were selected as feature variables. These, along with primary and secondary labels, resulted in a multivariate time series sample size of 144 * 27. The information within each multivariate time series is continuous in time, which is more conducive to uncovering short-term dependencies and cross-sectional correlations between variables. After data preprocessing, 8739 multivariate time series samples were obtained, including 7432 normal samples and 1307 fault samples obtained based on secondary labels. This indicates that these 1307 multivariate time series samples represent one or more types of faults. To reduce training pressure and alleviate data imbalance, 3000 samples were randomly selected from the normal samples and mixed with fault samples. These 4307 multivariate time series samples were then divided into training, validation, and test sets in a ratio of 0.7:0.15:0.15. Furthermore, the first layer of the stacked LSTM model had 100 hidden layers and 50 fully connected neurons, while the second layer had 50 hidden layers and 6 fully connected neurons. The training epochs were set to 200, and the batch size to 36. For diagnostic reliability and effectiveness validation, evaluation was first performed on a multivariate time series sample set basis, followed by an evaluation of the fault diagnosis accuracy at each time point within the multivariate time series, and finally, a comparison with existing methods was conducted.

[0118] (2) Evaluation of the effectiveness of the method:

[0119] First, the effectiveness and accuracy of the proposed method for fault classification are verified, including the sample aggregate diagnosis of multivariate time series samples based on secondary labels, as well as the fault diagnosis of samples at each time point based on primary labels. Figure 3 Fault classification results based on multivariate time series samples on the validation set and test set are presented respectively. It can be found that the stacked LSTM model shows good performance on both the validation set and the test set, with accuracies of 98.76% and 98.61% respectively.

[0120] Furthermore, to verify the effectiveness of the first-level label and probability correction method, two multivariate time series samples containing fault samples were randomly selected in the example. The states at each time step were diagnosed under the conditions of probability correction and without probability correction. Figure 4The confusion matrix is ​​used to show the state diagnosis results of two randomly selected multivariate time series samples before and after fault probability correction. It can be found that sample I contains three different fault labels, while sample II contains only two fault labels of the same class. After fault probability correction, the diagnosis results of the state at each time step within the two randomly selected samples are significantly improved, and the classification accuracy of sample II has reached 100%.

[0121] Considering the randomness of the above-mentioned random sampling, Table 2 also presents the average accuracy A based on the state diagnosis results at each time step within all multivariate time samples in the test set. sec And the average false positive rate index F NR Average false alarm rate index F PR Comparison of Indicators. As can be seen from the data listed in Table 2, the overall average accuracy is significantly improved after probability correction, and the corresponding average false negative rate and false positive rate are significantly reduced, thus demonstrating the effectiveness of the method proposed in this invention.

[0122] Table 2 Comparison of average evaluation metrics based on all test sets before and after probability correction.

[0123] Before probability correction 95.3％ 7.34％ 2.36％ After probability correction 98.2％ 3.27％ 1.14％

Claims

1. A wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration, characterized in that... Includes the following steps: Step 1: Analyze the characteristics of wind turbine SCADA data, preprocess the wind turbine SCADA data, and construct a multivariate time series sample set with secondary classification labels based on the matrix profile analysis method of time series. Step 2: Design a fault diagnosis model suitable for multivariate time series sample sets using a stacked LSTM model. This model performs fault diagnosis on the overall multivariate time series sample set and makes preliminary predictions on the sample state at each time step within the set. Step 3: Employ a hierarchical fault probability correction strategy to integrate multivariate information to correct the probability of samples belonging to different fault types at each time point, and select the class corresponding to the highest corrected probability value as the final fault diagnosis result. Step two includes the following steps: Step 2.1: Based on the stacked LSTM model, mine the dependencies between time series data and faults across multiple time periods. The input to the stacked LSTM model is a data matrix from the multivariate time series samples. X i As input, the first-level label vector Y i and secondary tags f i These represent the output training models; where the loss function corresponding to the first-level labels is the multivariate cross-entropy loss function. L m Second-level tags f i The corresponding loss function is the binary cross-entropy loss function. L b The corresponding expressions are shown in equations (5) and (6) respectively: （5）； （6）； In equations (5) and (6): The superscript indicates the number of multivariate time series samples; L represents the number of SCADA samples contained in each multivariate time series sample. This indicates that the corresponding category label is c; K Indicates the total number of tags contained within a first-level tag, subscript. Indicates the first i A multivariate time series sample, index t This represents the first multivariate time series sample collected by the SCADA system. t A sample at each moment; Indicate whether the label of the t-th SCADA sample in the i-th multivariate time series sample is c; This represents the probability that the model predicts the label of the t-th SCADA sample in the i-th multivariate time series sample to be c; Log() is the logarithm calculation; y i and p i These represent the categories when classified using second-level tags, specifically the first... i The labels of each sample are positive and their corresponding predicted probabilities; Step 2.2: Optimize the parameters of the stacked LSTM model using a grid search method and corresponding evaluation metrics; the evaluation metrics include accuracy metrics based on multivariate time series samples. A fir The average diagnostic accuracy of samples at each time point within all multivariate time series samples. A sec Average false positive rate F NR Average false alarm rate F PR ; （7）； （8）； （9）； （10）； In the formula: This represents the length of the sliding window used to construct multivariate time series samples, which includes... L Original samples collected by a SCADA system; N TP The number of normal samples identified as normal in a multivariate time series sample; N TN This represents the number of fault classes identified as fault classes in a multivariate time series sample. Indicates the first i The number of fault classes incorrectly identified as normal classes in all time-series samples within a multivariate time series sample; This represents the number of normal samples that are incorrectly identified as faulty samples within all time steps of the i-th multivariate time series sample.

2. The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration according to claim 1, characterized in that: Step one includes the following steps: Step 1.1: Preprocess the SCADA data of the wind turbine to eliminate outliers and missing values; Step 1.2: Combining the SCADA data and fault logs of the wind turbine, label the SCADA sampled data at each preprocessed moment; where no time is recorded in the fault log, it indicates normal operation and is labeled with the integer 0; other faults recorded in the fault log are labeled with different integers, and the corresponding sample label set is Y=[0,1,2…,K], where K is the number of fault categories; Step 1.3: Select strongly correlated monitoring quantities as feature variables, and use these selected strongly correlated monitoring quantities as the basic variables for constructing a multivariate time series sample set, as follows: Pearson correlation analysis is used to select characteristic variables, and the formula is as follows: （1）； In the formula: A and B represent the sets of candidate variables and the set of target variables, respectively; a i point b i Let i and b represent the i-th elements in A and B, respectively. and These are the means of the two feature subsets, respectively. P r This represents the Pearson correlation coefficient between two feature subsets; P r The value reflects the correlation between two feature subsets. Its range is [-1, 1]. A positive value indicates a positive correlation, a negative value indicates a negative correlation, and a value of 0 indicates no correlation. The closer the absolute value is to 1, the stronger the correlation between the two feature subsets. Step 1.4: Standardize the selected monitoring variables, including wind speed, wind direction, various rotational speeds, various temperature monitoring quantities, electrical monitoring quantities and output power, using the maximum and minimum value scaling shown in Equation (2); （2）； In formula (2): x For data to be processed, x * For the processed data, x max and x min These are the maximum and minimum values ​​of the monitored data for this variable, respectively. Step 1.5: Based on the selected feature variable data, a sliding window with a fixed window width is used to sample at a certain step size to obtain a multivariate time series data matrix containing multiple consecutive time points. X i Furthermore, with X i Construct label vectors for the corresponding labels at each time step. Y i As X i The first-level tag, and through the Y i Perform an OR operation on the tags inside to get X i The second-level labels, then the newly constructed multivariate time series sample set takes the following form: （3）； In formula (3): d Indicates the number of selected feature variables; L Indicates the length of the sampling window; N This represents the total number of multivariate time series samples; vector Y i The value of each element is one of [0, 1, 2, ..., K], where 0, 1, 2, ..., m represent the labels set according to the fault log, 0 represents no fault, and K represents the label of the Kth type of fault. f i This is a secondary label for multivariate fault samples. A value of 0 indicates that there is no fault during the sliding window sampling period, and a value of 1 indicates that there is a fault during the sliding window sampling period. This represents the i-th multivariate time series matrix sample constructed; Step 1.6: Divide the constructed multivariate time series into training set, validation set and test set according to a certain ratio. The training set is used for model training and parameter tuning, the validation set is used for performance evaluation of the trained model, and the test set contains several samples with known fault types for verifying the actual effect of the model.

3. The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration according to claim 2, characterized in that: In step 1.1, the 3δ criterion is used to initially identify outliers caused by sampling anomalies in various monitoring quantities, and interpolation is used to replace and supplement outliers and missing values ​​in the original data. The commonly used Lagrange interpolation method is used to handle outliers and missing values, as follows: Suppose we have n data points ( x 1, y 1),( x 2, y 2),…,( x n , y n Construct an n-1 degree polynomial P(x) at these data points, such that P( x i )= y i ( i =1,2,…,n); then the Lagrange interpolation formula for P(x) is: (4)； in, y i Indicates the first i The function value of each data point; x i Indicates the first i The independent variable values ​​for each data point; the π symbol in the formula represents the values ​​for all j ≠ i. j Perform product operations.

4. The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration according to claim 2, characterized in that: In step 1.5, the original SCADA system sampling data at each moment is called a SCADA sample. Each SCADA sample contains the values ​​of multiple monitoring quantities collected by different sensors at the same moment and corresponds to a primary label. Each multivariate time series sample contains L SCADA samples and their corresponding labels, as well as a secondary label of the multivariate time series sample.

5. The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration according to claim 1, characterized in that: In the indicators of step 2.2: A fir This indicates the classification accuracy of the stacked LSTM model for a multivariate time series sample population; A sec The average false positive rate (MFR) is the accuracy of the model in classifying samples across all time points within a multivariate time series sample. F NR This indicator reflects the missed detection rate of the model in fault diagnosis at each time step; the average false alarm rate metric. F PR This is used to reflect the false alarms in the model's fault diagnosis at each time point.

6. The wind turbine fault diagnosis method based on multivariate time series analysis and hierarchical probability calibration according to claim 1, characterized in that: Step three includes the following steps: Step 3.1: Calculate the global failure probability of multivariate time samples. P i As top-level prior information, the probability of individual failures The values ​​are corrected to eliminate the influence of individual extreme cases on the diagnostic results; an empirical weighting method is used to... The value is corrected once: （11）； In equation (11): For the first i Within a multivariate time series sample, the element at time t belongs to the... c The result of the first correction of the failure probability; For probability fusion weights, The value of reflects the confidence level of the model in predicting individual outcomes; To determine the expected values ​​of the probabilities of various types of faults under the condition of global sample faults in a multivariate time series; Indicates taking The superscript c with the highest probability value is the time c. t The maximum probability of fault classification in the sample remains consistent with the global probability of fault in the sample. After a correction using the top-level probability t The probability of fault classification for the sample at time step is , These represent the faults (0, 1, 2, ...) that the sampled at time t within the i-th multivariate time series sample belongs to after one correction. K The probability of ) is calculated; the corrected probability value takes into account the global information statistical distribution of multivariate time samples. Step 3.2: Fault probability estimation method based on statistical inference, including the following steps: Step 3.2.1: Based directly on the preprocessed SCADA sample data, construct separate systems according to category. K Different training sets D c c=0,1,2,… K The front of the dataset d The columns represent different monitoring variables, and the last column is the corresponding fault sample label; since the number of samples for each type may not be the same, D 0~ D K The number of rows in the matrix is ​​different, as shown in equation (12): D 0 and D c subscript l and m represent the number of normal type samples, respectively. l The first, and the second c There are m fault samples of this type; （12）； In equation (12): This represents the standardized value of the d-th monitoring variable in the m-th sample of the c-th type of fault; Step 3.2.2: A nonparametric kernel density estimation algorithm is used to learn from historical data of different fault types, to mine the conditional dependencies between various monitoring data and different fault labels, and to obtain the conditional probability density functions of various fault types. Then, given a new sample conditions, Let d represent the d different monitoring variables collected by the SCADA system at time t; using different We can estimate the probability that the sample belongs to different types of faults, and the probability is expressed as: （13）； In equation (13): Indicates at a given input time t Corresponding monitoring quantity In this case, utilize c Conditional probability density function of fault class It is estimated to belong to c The probability of a type of failure; m The number of training samples; h It is a hyperparameter, called bandwidth or window; d Indicates the number of monitored variables; The Gaussian kernel function; The training set is constructed from all sample data of the c-th type of fault. The meaning is the multiplication of d h's; This represents the input condition vector at a given time t. Indicated The condition vector acquired at the i-th time step; Estimating the probability density functions for different fault types under the same input conditions yields the probabilities belonging to different fault types. , Let represent the probabilities of a sample at time t belonging to the 0th, 1st, 2nd, ..., Kth time units, respectively, based on kernel density estimation. It estimates the probability using all samples of each type of fault in the training sample space. The estimation result contains the true distribution information of each type of fault sample and the individualized information of each sample. Therefore, this fault probability can be used as the underlying prior probability for further analysis. Make corrections; Step 3.3: Using Bayes' theorem, in Based on the conditional kernel density estimation results of various fault types As the underlying prior probability, for t The probability of time-based samples belonging to different faults is corrected a second time, specifically as follows: （14）； In equation (14): For the first i Within the multivariate time series samples, the one in the first t The time belongs to the first c The posterior probability of a fault class after secondary correction; The results obtained based on different kernel density estimates t The underlying prior probability that the sample at time step c belongs to fault c; when we obtain t The set of posterior probabilities that a sample at time step belongs to each class. back, Let and represent the probabilities that the sample at time t belongs to the 0th, 1st, 2nd, ..., Kth time groups after secondary correction; finally, the sample with the highest probability value is selected as the fault type to be diagnosed, and the corresponding reliability probability is . .

Images