Multivariable wind power curve modeling method based on SGBRT and GWO

Abstract

The present invention belongs to the technical field of power curve modeling, and specifically relates to a multivariable wind power curve modeling method based on SGBRT and GWO, the steps of which include dividing the operation area of ​​the wind turbine into five different stages, performing data cleaning according to the power control strategies of different stages; selecting input variables according to partial mutual information PMI, determining the correlation between external conditions and output power by comparing mutual information, and using PMI to process the coupling relationship between input variables; modeling based on stochastic gradient boosted regression tree SGBRT to obtain an SGBRT model; using the gray wolf optimization algorithm GWO to optimize and update the regression parameters in the SGBRT model; and modeling a multivariable wind power curve based on the optimized SGBRT model. The present invention realizes high-accuracy, high-efficiency and strong-adaptability wind turbine power curve prediction by combining SGBRT and GWO algorithms, and has significant performance advantages.

Description

Technical Field

[0001] The invention belongs to the technical field of power curve modeling, and in particular relates to a multivariable wind power curve modeling method based on SGBRT and GWO. Background Art

[0002] The power curve of a wind turbine (WT) can directly indicate its ability to effectively convert wind energy into electrical energy. Accurate modeling and monitoring of the power curve can achieve timely wind energy forecasting, early fault detection and condition monitoring of wind turbines, which can also be used to formulate effective control strategies and provide technical support for the smooth grid connection of large-scale wind power. Due to various practical factors such as location, wind speed, wind direction, ambient temperature, atmospheric pressure, machinery, control, deployment of wind farms, measurement errors and data processing methods, the empirical power curve may be significantly different from the theoretical power curve. Therefore, it is of great significance to use advanced algorithms to improve the modeling accuracy of wind turbine power curves, which is crucial to improve wind energy utilization and reduce overall operation and maintenance costs.

[0003] In existing research, the achievements of improved random nonparametric modeling methods in power curve models are remarkable, especially in terms of interpretability, accuracy and generalization ability. However, the ever-changing and complex operating environment of wind turbines introduces many external factors that may affect their power output, and ignoring certain input variables or their coupling effects may also reduce the accuracy of the power curve model. In addition, existing methods mainly focus on improving modeling accuracy, but as the modeling process of wavelet transform power curves develops towards the use of multivariate inputs, further research is still needed in terms of enhancing the adaptability of feature selection, training speed and robustness of the model to noisy data. Therefore, it is necessary to find a multivariate power curve modeling method based on advanced machine learning technology. Summary of the invention

[0004] In view of the above deficiencies in the prior art, the purpose of the present invention is to provide a multivariable wind power curve modeling method based on SGBRT and GWO. By combining the SGBRT and GWO algorithms, high-accuracy, high-efficiency and strong-adaptability wind turbine power curve prediction is achieved, which has significant performance advantages.

[0005] To achieve the above objectives, the technical solution adopted by the present invention is a multivariable wind power curve modeling method based on SGBRT and GWO, comprising the following steps:

[0006] S1. Divide the operation area of ​​the wind turbine into five different stages (i.e., divide the operation range of the wind turbine into five specific areas), including the cut-in stage, the maximum wind energy tracking stage, the transition stage, the rated power stage, and the cut-out stage, and perform data cleaning according to the power control strategies of different stages;

[0007] S2, selecting input variables according to partial mutual information PMI, determining the correlation between external conditions and output power by comparing mutual information, and using PMI to process the coupling relationship between input variables;

[0008] S3, modeling based on stochastic gradient boosted regression tree SGBRT to obtain the SGBRT model;

[0009] S4, using the Grey Wolf Optimization Algorithm GWO to optimize and update the regression parameters in the SGBRT model;

[0010] S5. Modeling of multivariable wind power curve based on the optimized SGBRT model.

[0011] As a preferred embodiment of the present invention, in S1, V is defined Wind is the actual wind speed; V Cut_in is the cut-in wind speed, i.e. the maximum operating wind speed of the fan during the cut-in stage; V Rated is the rated wind speed, i.e. the maximum operating wind speed of the wind turbine in the maximum wind energy tracking stage; V Cut_out is the cut-out wind speed, i.e. the maximum operating wind speed of the wind turbine at rated power stage; P Rated is the rated power of the wind turbine; a variable pitch control system is used to control and adjust the blade pitch angle and the wind rotor speed;

[0012] The various stages and their control strategies are:

[0013] Cut-in phase: When V Wind <V Cut_in When , the fan is turned off;

[0014] Maximum wind energy tracking stage: When V Cut_in ≤V Wind <V Rated When V Wind Continuously adjust to keep the fan operating at the optimal tip speed ratio;

[0015] Transition phase: When V Wind In V Rated When the wind speed fluctuates up and down, the pitch angle of the blade pitch angle is controlled between -5° and 5°, and the wind rotor speed is kept at the rated value to ensure the safety and stability of the wind turbine operation;

[0016] Rated power stage: When V Rated <V Wind ≤V Cut_out When the blade pitch angle is increased, the output power fluctuation is reduced;

[0017] Cut-out stage: When VWind >V Cut_out When the fan is turned off, the pitch angle is adjusted to the feathering state to stop the operation of the fan;

[0018] The data cleaning is as follows: the original data includes wind speed, output power, rotor speed, environmental conditions, transmission component temperature and vibration data. For the cut-in and cut-out stages, as the wind turbine stops working in the cut-in and cut-out stages, the non-zero output power value is cleared; for the maximum wind energy tracking stage, the transition stage and the rated power stage, the wind speed and output power scatter plots of these three stages show a band distribution. According to the corresponding power control strategy, the hybrid Archimedean-Copula function is used to establish an adaptive confidence boundary model, and the original data is cleaned in the three-dimensional Copula space.

[0019] As a preferred solution of the present invention, in S1, the specific process of data cleaning is:

[0020] S1.1. Define the wind speed in the maximum wind energy tracking stage, transition stage and rated power stage as V, and the output power as P, and estimate the cumulative marginal probability distribution of the wind speed and output power in the maximum wind energy tracking stage, transition stage and rated power stage (the estimation here can be done using existing technology), which are expressed as and , construct a joint frequency histogram;

[0021] S1.2. Fit the mixed Archimedean Copula function of each operating region, expressed as:

[0022] (1);

[0023] In the formula, u represents ;v stands for ; e is a natural constant; , , Represent Gumbel, Clayton and Frank Copula functions respectively; , , are the dependency parameters of the Gumbel, Clayton and FrankCopula functions, respectively, indicating the dependency relationship between random variables; Represents a Copula function; , , is a random number between 0 and 1; ln represents the logarithmic function with e as the base; exp represents the exponential function with e as the base;

[0024] S1.3, based on formula (1), derive exist The conditional probability under for:

[0025] (2);

[0026] In the formula, represents partial differential;

[0027] S1.4. The confidence level is set to , limiting the number of abnormal data, where is the complement of the confidence level; define an asymmetric coefficient α, and derive the confidence level through formula (3): The probability quantile of :

[0028] (3);

[0029] In the formula, , represent the low probability quantile and the high probability quantile respectively;

[0030] S1.5. Define a wind speed V0 and estimate its cumulative marginal probability distribution in the maximum wind energy tracking stage, transition stage and rated power stage. After that, Conditional The conditional probability becomes a single variable. Under wind speed V0, the output of formula (3) is , The corresponding , ; Then the quantile of P under V0 can be derived:

[0031] (4);

[0032] In the formula, Represents the inverse function of nonparametric conditional probability; , They represent the lower confidence limit and upper confidence limit of the wind power curve under V0 respectively;

[0033] S1.6. During the maximum wind energy tracking phase, transition phase, and rated power phase, the upper and lower limits of each wind speed are calculated at specified intervals and plotted in all operating areas based on these points to obtain the equivalent confidence boundaries of the wind power curve, through which abnormal data are identified and eliminated.

[0034] As a preferred solution of the present invention, in S2, the process of selecting input variables according to partial mutual information PMI is:

[0035] S2.1. For the random output vector y i , label the output vector set as , N is the number of samples, that is, the number of vectors; the information entropy of Y Calculated by formula (5):

[0036] (5);

[0037] In the formula, p i is corresponding to y i The probability distribution of

[0038] S2.2. For random input variable u i , marking its collection as ,u i Will affect y i , for the mutual information between U and Y, the elements contained in Y are calculated by the following formula:

[0039] (6);

[0040] In the formula, Indicates u i and i The joint probability density function of , Respectively represent u i ,y i The marginal probability density function of ; represents the mutual information between the two in the brackets (here it is the mutual information between U and Y, and the same applies to the rest); log is the logarithmic function, that is, the inverse function of the exponential function with base 10;

[0041] S2.3. For a sampling point u in U, its marginal probability density function (estimated density) It is expressed as:

[0042] (7);

[0043] Where d represents the dimension of the random input variable; is the sample covariance matrix, yes The determinant of ; for u and u i The Mahalanobis distance between them, h is the kernel density function, expressed as:

[0044] (8);

[0045] S2.4. If there is another input set U1, then when there is a mutual correlation between U and U1, The calculated value of cannot match the actual value. In this case, after removing the information about U contained in U1 and Y, PMI is needed to measure the correlation between U and Y. Once the information about U is eliminated, the initial U1 and Y are updated to U. T and Y T :

[0046] (9);

[0047] Where U T and Y T is the set of variables updated after removing the information about U contained in U1 and Y (i.e., the remaining information); represents the kernel regression estimate of U1 on U; represents the kernel regression estimate of Y on U, expressed as:

[0048] (10);

[0049] The PMI between S2.5, U1 and Y is calculated by the following formula:

[0050] (11);

[0051] In the formula, represents the partial mutual information between U1 and Y;

[0052] S2.6, define C as the total input variable set, X as the optimal input variable set, C X It represents the set of candidate input variables with the largest PMI value, and variable selection is performed based on PMI.

[0053] As a preferred solution of the present invention, in S2.6, the process of variable selection based on PMI is:

[0054] S2.6.1. Initially, , represents the empty set;

[0055] S2.6.2 If , construct the estimator according to formula (10) , represents the kernel regression estimate of Y on X;

[0056] S2.6.3. Calculate the remaining information ;

[0057] S2.6.4. For any element C in C j , according to formula (10), we can deduce , which means C j Kernel regression estimation on X;

[0058] S2.6.5. Calculate the remaining information ;

[0059] S2.6.6. Calculation , and select the largest C X ;

[0060] S2.6.7. Calculation and evaluation conditions:

[0061] (12);

[0062] Where P c represents the evaluation value; is the regression residual of Y calculated with the selected variables; p is the number of screened variables;

[0063] S2.6.8 If P c If the value decreases, C X Select into X;

[0064] S2.6.9, return to step S2.6.2 and continue the loop calculation. If P c The value no longer decreases, that is, P c When the value reaches its minimum, the selection of optimal input variables in X is completed.

[0065] As a preferred solution of the present invention, in S3, SGBRT is an improvement of the gradient boosting decision tree GBDT. When modeling, given a training sample set M, , x i is an element in the optimal input variable set, ; Calculate the average value of the output vector set Y for prediction of X, then calculate the residual of the current GBDT model, and fit a new decision tree to these residuals until the loss function has the minimum mathematical expectation value among all sample distributions, expressed as:

[0066] (13);

[0067] In the formula, Represents the mapping function from X to Y; represents the expected value of the joint distribution of Y and X; is the loss function, which represents the measure of the difference between the two in brackets (here Y and argmin represents the variable value that makes a function reach the minimum value.

[0068] As a preferred solution of the present invention, the algorithm flow of GBDT is:

[0069] S3.1. Define an initial constant function , minimize the loss function, and create a tree with the root node as the root node, expressed as:

[0070] (14);

[0071] Where β is the initial prediction value of the GBDT model;

[0072] S3.2, iteratively construct the K-boosting tree, including two loops. After the first loop "k=1 to K" and the second loop "i=1 to N", all regression trees are added together and the GBDT model is continuously updated, where K is the number of trees; the functional relationship calculated after iteration is as follows:

[0073] (15);

[0074] In the formula, is the expansion coefficient of the boosted tree; is the regression tree function; a k is the set of parameters for fitting the training data, which is expressed as ;

[0075] K-boosting tree refers to building a model by iteratively adding decision trees (enhanced trees). Each new tree attempts to correct the prediction error of the previous tree. The regression tree function is a decision tree that stores continuous values ​​instead of classification labels at the leaf nodes of the tree. This type of decision tree is used for regression problems, that is, predicting a continuous numerical output variable.

[0076] In order to solve the loss function in the regression tree function, the least squares method is used to fit the current pseudo residual, the negative gradient of the loss function The calculation is as follows:

[0077] (16);

[0078] In the formula, , Respectively represent from x i , mapping function from x to Y, where x is a sample point in X; represents the prediction function of the GBDT model for the element in the brackets before the kth iteration (here, it is the prediction function for the GBDT model x, and the same applies to the rest); k is the number of the current iteration, indicating that k decision trees have been added to the GBDT model; Used as an estimate of the residual;

[0079] Using a k and β k , jointly adapt to the training data, expressed as:

[0080] (17);

[0081] In the formula, is the regression tree function;

[0082] S3.3, using the least squares method to fit, a k It is expressed as:

[0083] (18);

[0084] Where a is a sample point in A;

[0085] β k The optimal value of is determined by:

[0086] (19);

[0087] S3.4. In each iteration from k=1 to K, the regression tree function used as a leaf node The input space is divided into L non-overlapping regions, represented as , represents the region corresponding to the lth leaf node of the kth tree, and the input space is the optimal input variable set X;

[0088] In each region, a single constant is predicted as follows:

[0089] (20);

[0090] Where I is an indicator function. When it is established, I=1, otherwise I is equal to 0; yes Estimated residuals in the region The mean function is a mathematical operation for finding the average value.

[0091] Therefore, in each sub-region In the regression tree, the lth leaf node is solved in the kth iteration. Since the regression tree obtained from equation (20) can predict the , so it is marked as The solution of equation (19) is simplified to:

[0092] (twenty one);

[0093] In the formula, Represents the predicted value of the lth leaf node (terminal node);

[0094] S3.5. In order to supervise the learning rate of the prediction process, all learners are normalized, and the learning rate adjustment factor z is set to adjust the step size of the gradient descent, 0<z≤1. The currently fitted GBDT model is updated in the following way:

[0095] (twenty two);

[0096] In the formula, Represents the prediction function of the GBDT model for X at the kth iteration;

[0097] S3.6. When the iteration ends, the final It is expressed as:

[0098] (twenty three).

[0099] As a preferred solution of the present invention, in S3, the improvement of SGBRT is that a randomization process is introduced into GBDT, that is, a random subset of the training data set is sampled to fit the base learner, and then the updated SGBRT model of the current iteration is calculated based on the separated random subset data set, and the steps include:

[0100] For a given training sample set M, the input space X is divided into L non-overlapping regions ;

[0101] definition is the loss function, is a set of randomly permuted integers The obtained sequence set, is the i-th sequence;

[0102] In each sample set of the disjoint region, a sample subset M is randomly selected without replacement s , denoted as , N s is the sample size and N s <N, , Corresponding to x i and i ;

[0103] Using M s Train the regression tree to obtain a mapping from X to Y; use the random subset as input to generate a new regression tree model for each region, and update equation (21) as:

[0104] (twenty four);

[0105] In the formula, represents the updated solution of equation (19);

[0106] The final SGBRT model is determined through iteration and is expressed as:

[0107] (25).

[0108] As a preferred solution of the present invention, in S4, the GWO algorithm is used to further optimize the regression parameters in the SGBRT model, including the number of trees K, the number of non-overlapping regions L and the learning rate adjustment factor z.

[0109] As a preferred solution of the present invention, in S4, the optimization process of the GWO algorithm is:

[0110] In S4.1, GWO algorithm, gray wolves are divided into four categories: I, II, III and IV. Wolf I is the chief leader, wolves II and III are subordinates, and the remaining IV wolves are responsible for the balance of relationships within the population. The position of the wolf pack is iteratively updated according to the fitness value of the wolf pack, and its updated position is used to search for the optimal solution in the problem space. During the hunting process, the behavior of gray wolves encircling prey is expressed as:

[0111] (26);

[0112] In the formula, t represents the current number of iterations; , Respectively represent the position vector of the gray wolf at the tth and t+1th iterations; represents the position vector of the prey at the tth iteration; Q and Z are coefficient vectors, respectively expressed as:

[0113] (27);

[0114] Where b1 and b2 are random constant vectors between 0 and 1; q is a convergence factor that decreases linearly from 2 to 0;

[0115] S4.2, when the prey position is identified, wolves II and III coordinate all wolves to surround the prey under the guidance of wolf I, which is expressed as:

[0116] (28);

[0117] In the formula, , , Respectively represent the position vectors of wolf I, wolf II and wolf III, expressed as:

[0118] (29);

[0119] In the formula, , , They represent the position vectors of wolf I, wolf II and wolf III at the tth iteration respectively; Q1, Q2, Q3 are the coefficient vectors of wolf I, wolf II and wolf III respectively; , , Respectively represent the distances between wolf I, wolf II, and wolf III and other individuals, expressed as:

[0120] (30);

[0121] Where Z1, Z2, and Z3 are random vectors;

[0122] S4.3. Constructing the fitness function of GWO , expressed as:

[0123] (31);

[0124] In the formula, RMSE is the root mean square error; MAE is the mean absolute error; R 2 is the coefficient of determination; , , They are RMSE, MAE and R 2 A weighting factor for ; c is a positive constant used to prevent division by zero errors.

[0125] The algorithm involved in the present invention can be executed by an electronic device, which includes a memory, a processor, and a computer program stored in the memory and executable on the processor, and the above algorithm calculation is implemented by executing the software by the processor.

[0126] The beneficial effects of the present invention are:

[0127] This paper proposes a wind turbine power curve prediction model that combines the stochastic gradient boosting regression tree (SGBRT) and the grey wolf optimization (GWO) algorithm. By using the probability description method in the three-dimensional Copula space to clean up the original data, invalid and unnatural data are effectively eliminated. At the same time, partial mutual information analysis is used to determine the key input variables, which improves the data quality and model accuracy. The model shows excellent accuracy and efficiency in predicting the power curve of variable speed constant frequency wind turbines.

[0128] The proposed SGBRT-GWO model architecture performs well in terms of modeling accuracy, efficiency, and adaptability. Compared with parameterized, artificial neural network, and XGBoost models, the proposed SGBRT-GWO model architecture performs well in terms of normalized residual error, RMSE, MAE, and R 2The model iteration process is also simplified, reducing the computational cost and time while maintaining high prediction accuracy. These advantages make the scheme have important practical application value in optimizing wind energy production and improving the efficiency of wind turbines, and provide strong technical support for the operation, maintenance and control of wind power generation systems. BRIEF DESCRIPTION OF THE DRAWINGS

[0129] Figure 1 It is a schematic diagram of the process of the present invention;

[0130] Figure 2 It is a schematic diagram of the variable selection process of the PMI algorithm in the verification process of the present invention;

[0131] Figure 3 It is a schematic diagram of the confidence boundary in the verification process of the present invention;

[0132] Figure 4 It is a schematic diagram of the results of the SGBRT-GWO prediction power curve model during the verification process of the present invention;

[0133] Figure 5 It is a schematic diagram of the normalized prediction residual of the prediction power curve model in the verification process of the present invention;

[0134] Figure 6 It is a schematic diagram of RMSE in each iteration during the verification process of the present invention;

[0135] Figure 7 yes Figure 6 A magnified view of the structure at center A;

[0136] Figure 8 It is a schematic diagram of the prediction results of the five-parameter method in the verification process of the present invention;

[0137] Fig. 9 It is a schematic diagram of the prediction results of the improved ANN in the verification process of the present invention;

[0138] Fig.10 It is a schematic diagram of the prediction results of the XGBoost algorithm in the verification process of the present invention;

[0139] Fig.11 It is a schematic diagram of the prediction results of SGBRT-GWO in the verification process of the present invention;

[0140] Fig.12 It is a schematic diagram of MAE comparison of different power curve models in the verification process of the present invention;

[0141] Fig.13 is the R of different power curve models in the verification process of the present invention 2 Comparison diagram. DETAILED DESCRIPTION

[0142] The embodiments of the present invention are further described below in conjunction with the accompanying drawings:

[0143] like Figure 1 As shown, the multivariable wind power curve modeling method based on SGBRT and GWO includes the following steps:

[0144] S1. Divide the operation area of ​​the wind turbine into five different stages, including the cut-in stage, the maximum wind energy tracking stage, the transition stage, the rated power stage and the cut-out stage, and perform data cleaning according to the power control strategies of different stages;

[0145] S2, selecting input variables according to partial mutual information PMI, determining the correlation between external conditions and output power by comparing mutual information, and using PMI to process the coupling relationship between input variables;

[0146] S3, modeling based on stochastic gradient boosted regression tree SGBRT to obtain the SGBRT model;

[0147] S4, using the Grey Wolf Optimization Algorithm GWO to optimize and update the regression parameters in the SGBRT model;

[0148] S5. Modeling of multivariable wind power curve based on the optimized SGBRT model.

[0149] In S1, define V Wind is the actual wind speed; V Cut_in is the cut-in wind speed, i.e. the maximum operating wind speed of the fan during the cut-in stage; V Rated is the rated wind speed, i.e. the maximum operating wind speed of the wind turbine in the maximum wind energy tracking stage; V Cut_out is the cut-out wind speed, i.e. the maximum operating wind speed of the wind turbine at rated power stage; P Rated is the rated power of the wind turbine; a variable pitch control system is used to control and adjust the blade pitch angle and the wind rotor speed;

[0150] The various stages and their control strategies are:

[0151] Cut-in phase: When V Wind <V Cut_in When , the fan is turned off;

[0152] Maximum wind energy tracking stage: When V Cut_in ≤V Wind <V Rated When V Wind Continuously adjust to keep the fan operating at the optimal tip speed ratio;

[0153] Transition phase: When V Wind In V RatedWhen the wind speed fluctuates up and down, the pitch angle of the blade pitch angle is controlled between -5° and 5°, and the wind rotor speed is kept at the rated value to ensure the safety and stability of the wind turbine operation;

[0154] Rated power stage: When V Rated <V Wind ≤V Cut_out When the blade pitch angle is increased, the output power fluctuation is reduced;

[0155] Cut-out stage: When V Wind >V Cut_out When the fan is turned off, the pitch angle is adjusted to the feathering state to stop the operation of the fan;

[0156] The division of the five stages can also be made with reference to other existing technologies. The pitch control system is a system for adjusting the angle of the blades of a wind turbine generator set to optimize the wind energy conversion efficiency and protect the wind turbine from extreme wind speeds. It can be implemented using known technologies.

[0157] The data cleaning is as follows: the original data includes wind speed, output power, rotor speed, environmental conditions, transmission component temperature and vibration data. For the cut-in and cut-out stages, as the wind turbine stops working in the cut-in and cut-out stages, the non-zero output power value is cleared; for the maximum wind energy tracking stage, the transition stage and the rated power stage, the wind speed and output power scatter plots of these three stages show a band distribution. According to the corresponding power control strategy, the hybrid Archimedean-Copula function is used to establish an adaptive confidence boundary model, and the original data is cleaned in the three-dimensional Copula space.

[0158] In S1, the specific process of data cleaning is as follows:

[0159] S1.1. Define the wind speed in the maximum wind energy tracking stage, transition stage and rated power stage as V, and the output power as P. Estimate the cumulative marginal probability distribution of wind speed and output power in the maximum wind energy tracking stage, transition stage and rated power stage, which are expressed as and , construct a joint frequency histogram;

[0160] S1.2. Fit the mixed Archimedean Copula function of each operating region, expressed as:

[0161] (1);

[0162] In the formula, u represents ;v stands for ; e is a natural constant; , , Represent Gumbel, Clayton and Frank Copula functions respectively; , , are the dependency parameters of the Gumbel, Clayton and FrankCopula functions, respectively, indicating the dependency relationship between random variables; Represents a Copula function; , , is a random number between 0 and 1;

[0163] The known parameters in formula (1) can be estimated using the expectation maximization method with gradient descent optimization;

[0164] S1.3, based on formula (1), derive exist The conditional probability under for:

[0165] (2);

[0166] In the formula, represents partial differential; this embodiment uses non-parametric conditional probability The method solves the quantile of P under a specific V;

[0167] S1.4. The confidence level is set to , limiting the number of abnormal data, where is the complement of the confidence level; define an asymmetric coefficient α, and derive the confidence level through formula (3): The probability quantile of :

[0168] (3);

[0169] In the formula, , denote low and high probability quantiles (lower and higher probability quantiles), respectively;

[0170] S1.5. Define a wind speed V0 and estimate its cumulative marginal probability distribution in the maximum wind energy tracking stage, transition stage and rated power stage. After that, Conditional The conditional probability becomes a single variable. Under wind speed V0, the output of formula (3) is , The corresponding , ; Then the quantile of P under V0 can be derived:

[0171] (4);

[0172] In the formula, Represents the inverse function of nonparametric conditional probability; , They represent the lower confidence limit and upper confidence limit of the wind power curve under V0 respectively;

[0173] S1.6. During the maximum wind energy tracking phase, transition phase, and rated power phase, the upper and lower limits of each wind speed are calculated at specified intervals and plotted in all operating areas based on these points to obtain the equivalent confidence boundaries of the wind power curve, through which abnormal data are identified and eliminated.

[0174] and The Copula conditional probability connection derived from Equation (4) can obtain the expected confidence level by comprehensively considering the joint and local probability distributions.

[0175] Since many external factors can affect the output power of WT (wind turbine / wind generator / wind turbine for short), it is crucial to choose appropriate input variables for power curve modeling, which is conducive to balancing the calculation dimension and accuracy.

[0176] In S2, the process of selecting input variables based on partial mutual information PMI is:

[0177] S2.1. For the random output vector y i , label the output vector set as , N is the number of samples, that is, the number of vectors; the information entropy of Y Calculated by formula (5):

[0178] (5);

[0179] In the formula, p i is corresponding to y i The probability distribution of

[0180] S2.2. For random input variable u i , marking its collection as ,u i Will affect y i , for the mutual information between U and Y, the elements contained in Y are calculated by the following formula:

[0181] (6);

[0182] In the formula, Indicates u i and i The joint probability density function of , Respectively represent u i ,y i The marginal probability density function of ; represents the mutual information between the two in brackets;

[0183] S2.3. For a sampling point u in U, its marginal probability density function It is expressed as:

[0184] (7);

[0185] Where d represents the dimension of the random input variable; is the sample covariance matrix, yes The determinant of ; for u and u i The Mahalanobis distance between them, h is the kernel density function, expressed as:

[0186] (8);

[0187] S2.4. If there is another input set U1, then when there is a mutual correlation between U and U1, The calculated value of cannot match the actual value. In this case, after removing the information about U contained in U1 and Y, PMI is needed to measure the correlation between U and Y. Once the information about U is eliminated, the initial U1 and Y are updated to U. T and Y T :

[0188] (9);

[0189] Where U T and Y T is the set of variables updated after removing the information about U contained in U1 and Y; represents the kernel regression estimate of U1 on U; represents the kernel regression estimate of Y on U, expressed as:

[0190] (10);

[0191] The PMI between S2.5, U1 and Y is calculated by the following formula:

[0192] (11);

[0193] In the formula, represents the partial mutual information between U1 and Y;

[0194] S2.6, define C as the total input variable set, X as the optimal input variable set, C X It represents the set of candidate input variables with the largest PMI value, and variable selection is performed based on PMI.

[0195] In S2.6, the process of variable selection based on PMI is:

[0196] S2.6.1. Initially, , represents the empty set;

[0197] S2.6.2 If , construct the estimator according to formula (10) , represents the kernel regression estimate of Y on X;

[0198] S2.6.3. Calculate the remaining information ;

[0199] S2.6.4. For any element C in C j , according to formula (10), we can deduce , which means C j Kernel regression estimation on X;

[0200] S2.6.5. Calculate the remaining information ;

[0201] S2.6.6. Calculation , and select the largest C X ;

[0202] S2.6.7. Calculation and evaluation conditions:

[0203] (12);

[0204] Where P c represents the evaluation value; is the regression residual of Y calculated with the selected variables; p is the number of screened variables;

[0205] S2.6.8 If P c If the value decreases, C X Select into X;

[0206] S2.6.9, return to step S2.6.2 and continue the loop calculation. If P c The value no longer decreases, that is, P c When the value reaches its minimum, the selection of optimal input variables in X is completed.

[0207] In S3, SGBRT is an improvement of the gradient boosted decision tree GBDT. When modeling, given a training sample set M, , x i is an element in the optimal input variable set, ; Calculate the average value of the output vector set Y for prediction of X, then calculate the residual of the current GBDT model, and fit a new decision tree to these residuals until the loss function has the minimum mathematical expectation value among all sample distributions, expressed as:

[0208] (13);

[0209] In the formula, Represents the mapping function from X to Y; represents the expected value of the joint distribution of Y and X; is the loss function, which represents the measure of the difference between the two in brackets.

[0210] The algorithm flow of GBDT is:

[0211] S3.1. Define an initial constant function , minimize the loss function, and create a tree with the root node as the root node, expressed as:

[0212] (14);

[0213] Where β is the initial prediction value of the GBDT model;

[0214] β is a parameter that needs to be optimized. It is a parameter that needs to be determined by minimizing the loss. In the first step It is defined as a constant function, and this constant value is β obtained by minimizing the loss function; this initial value (constant value) is the basis for subsequent iterations. In each subsequent step, a new decision tree will be constructed to correct this initial value.

[0215] S3.2, iteratively construct the K-boosting tree, including two loops. After the first loop "k=1 to K" and the second loop "i=1 to N", all regression trees are added together and the GBDT model is continuously updated, where K is the number of trees; the functional relationship calculated after iteration is as follows:

[0216] (15);

[0217] In the formula, is the expansion coefficient of the boosted tree; is the regression tree function; a k is the set of parameters for fitting the training data, which is expressed as ;

[0218] In order to solve the loss function in the regression tree function, the least squares method is used to fit the current pseudo residual, the negative gradient of the loss function The calculation is as follows:

[0219] (16);

[0220] In the formula, , Respectively represent from x i , mapping function from x to Y, where x is a sample point in X; represents the prediction function of the GBDT model for the elements in brackets before the kth iteration; k is the number of the current iteration, indicating that k decision trees have been added to the GBDT model; Used as an estimate of the residual;

[0221] Using a k and β k , jointly adapt to the training data, expressed as:

[0222] (17);

[0223] In the formula, is the regression tree function;

[0224] S3.3, using the least squares method to fit, a k It is expressed as:

[0225] (18);

[0226] Where a is a sample point in A;

[0227] β k The optimal value of is determined by:

[0228] (19);

[0229] S3.4. In each iteration from k=1 to K, the regression tree function used as a leaf node The input space is divided into L non-overlapping regions, represented as , represents the region corresponding to the lth leaf node of the kth tree, and the input space is the optimal input variable set X;

[0230] In each region, a single constant is predicted as follows:

[0231] (20);

[0232] Where I is an indicator function. When it is established, I=1, otherwise I is equal to 0; yes Estimated residuals in the region The average value of

[0233] Therefore, in each sub-region In the regression tree, the lth leaf node is solved in the kth iteration. Since the regression tree obtained from equation (20) can predict the , so it is marked as The solution of equation (19) is simplified to:

[0234] (twenty one);

[0235] In the formula, Represents the predicted value of the lth leaf node;

[0236] S3.5. In order to supervise the learning rate of the prediction process, all learners are normalized, and the learning rate adjustment factor z is set to adjust the step size of the gradient descent, 0<z≤1. The currently fitted GBDT model is updated in the following way:

[0237] (twenty two);

[0238] In the formula, Represents the prediction function of the GBDT model for X at the kth iteration;

[0239] S3.6. When the iteration ends, the final It is expressed as:

[0240] (twenty three).

[0241] In S3, the improvement of SGBRT is to introduce a randomization process in GBDT, that is, sampling a random subset of the training data set to fit the base learner, and then calculating the updated SGBRT model of the current iteration based on the separated random subset of the data set. This improvement reduces the possibility of overfitting and speeds up the algorithm by reducing the computational burden. The steps include:

[0242] For a given training sample set M, the input space X is divided into L non-overlapping regions ;

[0243] definition is the loss function, is a set of randomly permuted integers The obtained sequence set, is the i-th sequence;

[0244] In each sample set of the disjoint region, a sample subset M is randomly selected without replacements , denoted as , N s is the sample size and N s <N, , Corresponding to x i and i ;

[0245] Using M s Train the regression tree to obtain a mapping from X to Y; use the random subset as input to generate a new regression tree model for each region, and update equation (21) as:

[0246] (twenty four);

[0247] In the formula, represents the updated solution of equation (19);

[0248] The final SGBRT model is determined through iteration and is expressed as:

[0249] (25).

[0250] In S4, the GWO algorithm is used to further optimize the regression parameters in the SGBRT model, including the number of trees K, the number of disjoint regions L, and the learning rate adjustment factor z.

[0251] In the GWO algorithm, gray wolves are divided into four categories: I, II, III, and IV. Wolf I is the chief leader, wolves II and III are subordinates, and the remaining IV wolves are responsible for the balance of relationships within the population. The position of the wolf pack is iteratively updated according to the fitness value of the wolf pack, and its updated position is used to search for the optimal solution in the problem space.

[0252] The verification process is:

[0253] The practicality and advancement of the power curve model proposed in this embodiment are verified by using the SCADA measured data of a 1.5MW doubly-fed wind turbine in a wind farm. The time span of the data is four months. A total of 17,537 data groups are measured, recorded every 10 minutes, and each group contains 29 variables.

[0254] According to the PMI algorithm in this embodiment, variable selection is performed, and the schematic diagram of the variable selection process is shown in FIG. Figure 2 As shown. c The value decreases as the number of input variables increases. When 8 input variables are screened, P c The value reaches the minimum value (-1.805×10 4). These variables are selected based on the maximum PMI value in each screening process. The variables finally selected are wind speed, wind direction, pitch angle, yaw error, rotor speed, tip speed ratio, gearbox oil temperature, and generator temperature. Among the selected variables, the tip speed ratio is calculated using data such as wind turbine radius, rotor speed, and wind speed, and is used for PMI selection.

[0255] When cleaning data, set the confidence level to =0.95, according to the data cleaning method in the embodiment, determine the upper and lower confidence limits of the output power at different wind speeds , Using SCADA data, a schematic diagram of the confidence boundary is shown below: Figure 3 shown.

[0256] Through the data cleaning strategy proposed in this embodiment, non-zero output power values ​​are eliminated in the cut-in and cut-out stages. In the maximum wind energy tracking stage, transition stage, and rated power stage, by setting the upper and lower confidence boundaries of the power curve, abnormal data can also be efficiently identified and removed. According to the probability distribution analysis, compared with the original data set, the cleaned data group does not have a peak near 0kW power. Relatively speaking, the cleaned data is more concentrated in the normal power range, which is conducive to ensuring the accuracy and efficiency of subsequent power curve prediction.

[0257] After data cleaning, 11761 data sets are retained. To prevent underfitting, the training samples should account for 2 / 3-4 / 5 of the total data. In this case, the first 8761 sets of data are used to train the model, and the remaining 3000 sets of data are used as validation and test samples. The initial regression parameters K, L and z are set to 174, 47 and 0.05 respectively, and N is set to s / N=0.5, the sample data is input into the SGBRT-GWO model framework constructed in this embodiment, and the predicted power curve model (i.e., the SGBRT-GWO model framework) verified by the last 3000 data sets is trained. The results are as follows Figure 4 As shown, the blue star points represent the actual measurement data, and the red points are the power predicted by the model based on this embodiment. Figure 4 This indicates that the prediction model is not a single curve but a region, which means that for the same wind speed, there may be different predicted power output values. The reason is that the model proposed in this embodiment incorporates multiple input variables, making it closer to the actual output power characteristics of the WT.

[0258] In the regression process of the predicted power curve model, the closer the correlation coefficient value is to 1, the higher the fit is. During the validation process, the correlation coefficients of the training, validation, and test data were 0.99638, 0.99640, and 0.99617, respectively. The correlation coefficient between the predicted power curve model and the entire data set was 0.99635, confirming the excellent regression performance.

[0259] In order to quantify the modeling error, the normalized prediction residual of the predicted power curve model was calculated, as Figure 5 As shown, the RMSE in each iteration is also calculated, using the last 3000 sets of data as samples, as shown in Figure 6 and Figure 7 shown.

[0260] from Figure 5-Figure 7 It can be seen that under various external conditions, the normalized residual between the predicted power value and the actual power value is very small (within ±6.56%). As the decision tree iteration proceeds, the RMSE of the predicted power model for the training, validation, and test data sets gradually converges. After the 174th iteration, the RMSE stabilizes at about 27kW and remains consistent. In addition, the calculated final MAE is also very small, only 19.27 kW, while the R of the tested 1.5 MW wind turbine is 2 The value is also greater than 99.34%. These preliminary verifications show that the multivariable model proposed in this embodiment can effectively predict the power characteristics of wind turbines with satisfactory comprehensive accuracy.

[0261] To further verify the superiority of the proposed method, three existing modeling methods, including the five-parameter method, the improved ANN, and the XGBoost algorithm, were also used to model the power curve of the tested WT. The same clean validation and test data were used as input for all models. The comparison results are shown in Figure 2. Figure 8-Figure 11 As shown. The prediction model of the five-parameter method only uses wind speed as the system input, resulting in the predicted power following a smooth curve and failing to reflect the dynamic WT power characteristics with other external influencing factors. ANN, XGBoost, and the SGBRT-GWO-based model proposed in this embodiment all consider multiple input variables, resulting in a band-like distribution of the power curve.

[0262] A comparative analysis of the normalized prediction residuals of each model was performed, and the maximum residuals of the power curve model developed using the five-parameter method, improved ANN, XGBoost, and the proposed SGBRT-GWO were 25.99%, 15.32%, 14.44%, and 6.56%, respectively. The results show that the SGBRT-GWO model proposed in this embodiment is superior to other models in terms of residuals, which is mainly due to its comprehensive consideration of external factors and the optimization effect of GWO on regression parameters. The final RMSE of the established SGBRT-GWO model is only 27kW, which is about 25%, 45%, and 50% of the parameter-based, ANN, and XGBoost models. Compared with other non-parametric methods, the iteration process can also be reduced (174 trees) and the number of iterations exceeds 700 times.

[0263] Since some parts of the power curve model may have abnormal predictions, in order to ensure the overall accuracy of the predicted power curve, this verification example also calculates and compares the MAE and R of different models in different wind speed ranges (from 3m / s to 16m / s, with an interval of 0.5m / s). 2 ,like Fig.12 and Fig.13 The comparison results show that the MAE and R 2 The accuracy is improved. In different wind speed areas, the maximum average absolute power of the SGBRT-GWO model is 27.92kW, which is only 1.86% of the rated power (1.5MW), while the maximum average absolute powers of the models based on parameters, ANN, and XGBoost are 101.24kW (6.75%), 70.17kW (4.68%), and 63.45kW (4.23%), respectively. The proposed model has a good performance in most wind speed areas. 2 The maximum value can reach 98%, and the average is 98.61%. In comparison, the R 2 The value is much smaller and fluctuates greatly in different speed ranges. Average R of the models based on parameters, ANN and XGBoost 2 The values ​​are 71.91%, 92.45% and 92.49% respectively.

[0264] In summary, this embodiment uses a multivariate approach, combined with SGBRT and GWO algorithms, to develop an advanced power curve prediction model for variable speed constant frequency wind turbines. The proposed model is verified in actual operation using real SCADA data from WT. The comparison results verify the effectiveness and superiority of the proposed method. The specific conclusions include:

[0265] Through data cleaning of the confidence bound model, the original data is more concentrated in the normal operating area, and using the input variable selection of PMI, an ideal balance between calculation dimension and prediction accuracy can also be achieved.

[0266] The SGBRT-GWO based model shows impressive accuracy and efficiency in predicting the power curve of variable speed constant frequency wind turbines. The normalized residual error is less than 6.56% and the RMSE is about 27kW. In addition, the MAE is also quite low, accounting for only 1.86% of the rated power, R 2 The values ​​are also promising, reaching 98% in most wind speed areas and an average of 98.61%.

[0267] The improvement of the SGBRT-GWO based method is also verified. Compared with the parameter-based, ANN-based, and XGBoost-based models, the proposed power curve model performs better in all normalized residuals, RMSE, MAE, and R 2 In addition, the iteration process is also simplified.

Claims

1. A multivariable wind power curve modeling method based on SGBRT and GWO, characterized by The following steps are involved: S1. Divide the operation area of ​​the wind turbine into five different stages, including the cut-in stage, the maximum wind energy tracking stage, the transition stage, the rated power stage and the cut-out stage, and perform data cleaning according to the power control strategies of different stages; S2, selecting input variables according to partial mutual information PMI, determining the correlation between external conditions and output power by comparing mutual information, and using PMI to process the coupling relationship between input variables; S3, modeling based on stochastic gradient boosted regression tree SGBRT to obtain the SGBRT model; S4, using the Grey Wolf Optimization Algorithm GWO to optimize and update the regression parameters in the SGBRT model; S5. Modeling of multivariable wind power curve based on the optimized SGBRT model; In S2, the process of selecting input variables according to partial mutual information PMI is: S2.

1. For the random output vector y i , label the output vector set as , N is the number of samples, that is, the number of vectors; the information entropy of Y Calculated by formula (5): （5）； In the formula, p i is corresponding to y i The probability distribution of S2.

2. For random input variable u i , marking its collection as ,u i Will affect y i , for the mutual information between U and Y, the elements contained in Y are calculated by the following formula: （6）； In the formula, Indicates u i and i The joint probability density function of , Respectively represent u i ,y i The marginal probability density function of ; represents the mutual information between the two in brackets; S2.

3. For a sampling point u in U, its marginal probability density function It is expressed as: （7）； Where d represents the dimension of the random input variable; is the sample covariance matrix, yes The determinant of ; for u and u i The Mahalanobis distance between them, h is the kernel density function, expressed as: （8）； S2.

4. If there is another input set U1, then when there is a mutual correlation between U and U1, The calculated value of cannot match the actual value. In this case, PMI is needed to measure the correlation between U and Y after removing the information about U contained in U1 and Y. Once the information of U is eliminated, the initial U1 and Y are updated to U T and Y T : （9）； Where U T and Y T is the set of variables updated after removing the information about U contained in U1 and Y; represents the kernel regression estimate of U1 on U; represents the kernel regression estimate of Y on U, expressed as: （10）； The PMI between S2.5, U1 and Y is calculated by the following formula: （11）； In the formula, represents the partial mutual information between U1 and Y; S2.6, define C as the total input variable set, X as the optimal input variable set, C X It represents the set of candidate input variables with the largest PMI value, and variable selection is performed based on PMI.

2. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 1 is characterized in that: In S1, V is defined Wind is the actual wind speed; V Cut_in is the cut-in wind speed, i.e. the maximum operating wind speed of the fan during the cut-in stage; V Rated is the rated wind speed, i.e. the maximum operating wind speed of the wind turbine in the maximum wind energy tracking stage; V Cut_out is the cut-out wind speed, i.e. the maximum operating wind speed of the wind turbine at rated power stage; P Rated is the rated power of the wind turbine; a variable pitch control system is used to control and adjust the blade pitch angle and the wind rotor speed; The various stages and their control strategies are: Cut-in phase: When V Wind <V Cut_in When , the fan is turned off; Maximum wind energy tracking stage: When V Cut_in ≤V Wind <V Rated When V Wind Continuously adjust to keep the fan operating at the optimal tip speed ratio; Transition phase: When V Wind In V Rated When the wind speed fluctuates up and down, the pitch angle of the blade pitch angle is controlled between -5° and 5°, and the wind rotor speed is kept at the rated value to ensure the safety and stability of the wind turbine operation; Rated power stage: When V Rated <V Wind ≤V Cut_out When the blade pitch angle is increased, the output power fluctuation is reduced; Cut-out stage: When V Wind >V Cut_out When the fan is turned off, the pitch angle is adjusted to the feathering state to stop the operation of the fan; The data cleaning is as follows: the original data includes wind speed, output power, rotor speed, environmental conditions, transmission component temperature and vibration data. For the cut-in and cut-out stages, as the wind turbine stops working in the cut-in and cut-out stages, the non-zero output power value is cleared; for the maximum wind energy tracking stage, the transition stage and the rated power stage, the wind speed and output power scatter plots of these three stages show a band distribution. According to the corresponding power control strategy, the hybrid Archimedean-Copula function is used to establish an adaptive confidence boundary model, and the original data is cleaned in the three-dimensional Copula space.

3. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 2 is characterized in that: In the above S1, the specific process of data cleaning is as follows: S1.

1. Define the wind speed in the maximum wind energy tracking stage, transition stage and rated power stage as V, and the output power as P. Estimate the cumulative marginal probability distribution of wind speed and output power in the maximum wind energy tracking stage, transition stage and rated power stage, which are expressed as and , construct a joint frequency histogram; S1.

2. Fit the mixed Archimedean Copula function of each operating region, expressed as: （1）； In the formula, u represents ;v stands for ; e is a natural constant; , , Represent Gumbel, Clayton and Frank Copula functions respectively; , , are the dependency parameters of the Gumbel, Clayton and FrankCopula functions, respectively, indicating the dependency relationship between random variables; Represents a Copula function; , , is a random number between 0 and 1; S1.3, based on formula (1), derive exist The conditional probability under for: （2）； In the formula, represents partial differential; S1.

4. The confidence level is set to , limiting the number of abnormal data, where is the complement of the confidence level; define an asymmetric coefficient α, and derive the confidence level through formula (3): The probability quantile of : （3）； In the formula, , represent the low probability quantile and the high probability quantile respectively; S1.

5. Define a wind speed V0 and estimate its cumulative marginal probability distribution in the maximum wind energy tracking stage, transition stage and rated power stage. After that, Conditional The conditional probability becomes a single variable. Under wind speed V0, the output of formula (3) is , The corresponding , ; Then the quantile of P under V0 can be derived: （4）； In the formula, Represents the inverse function of nonparametric conditional probability; , They represent the lower confidence limit and upper confidence limit of the wind power curve under V0 respectively; S1.

6. During the maximum wind energy tracking phase, transition phase, and rated power phase, the upper and lower limits of each wind speed are calculated at specified intervals and plotted in all operating areas based on these points to obtain the equivalent confidence boundaries of the wind power curve, through which abnormal data are identified and eliminated.

4. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 1 is characterized in that: In S2.6, the process of variable selection based on PMI is: S2.6.

1. Initially, , represents the empty set; S2.6.2 If , construct the estimator according to formula (10) , represents the kernel regression estimate of Y on X; S2.6.

3. Calculate the remaining information ; S2.6.

4. For any element C in C j , according to formula (10), we can deduce , which means C j Kernel regression estimation on X; S2.6.

5. Calculate the remaining information ; S2.6.

6. Calculation , and select the largest C X ; S2.6.

7. Calculation and evaluation conditions: （12）； Where P c represents the evaluation value; is the regression residual of Y calculated with the selected variables; p is the number of screened variables; S2.6.8 If P c If the value decreases, C X Select into X; S2.6.9, return to step S2.6.2 and continue the loop calculation. If P c The value no longer decreases, that is, P c When the value reaches its minimum, the selection of optimal input variables in X is completed.

5. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 1 is characterized in that: In S3, SGBRT is an improvement of the gradient boosted decision tree GBDT. When modeling, given a training sample set M, , x i is an element in the optimal input variable set, ; Calculate the average value of the output vector set Y for prediction of X, then calculate the residual of the current GBDT model, and fit a new decision tree to these residuals until the loss function has the minimum mathematical expectation value among all sample distributions, expressed as: （13）； In the formula, Represents the mapping function from X to Y; represents the expected value of the joint distribution of Y and X; is the loss function, which represents the measure of the difference between the two in brackets.

6. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 5 is characterized in that: The algorithm flow of GBDT is: S3.

1. Define an initial constant function , minimize the loss function, and create a tree with the root node as the root node, expressed as: （14）； Where β is the initial prediction value of the GBDT model; S3.2, iteratively construct the K-boosting tree, including two loops. After the first loop "k=1 to K" and the second loop "i=1 to N", all regression trees are added together and the GBDT model is continuously updated, where K is the number of trees; the functional relationship calculated after iteration is as follows: （15）； In the formula, is the expansion coefficient of the boosted tree; is the regression tree function; a k is the set of parameters for fitting the training data, which is expressed as ; In order to solve the loss function in the regression tree function, the least squares method is used to fit the current pseudo residual, the negative gradient of the loss function The calculation is as follows: （16）； In the formula, , Respectively represent from x i , mapping function from x to Y, where x is a sample point in X; represents the prediction function of the GBDT model for the elements in brackets before the kth iteration; k is the number of the current iteration, indicating that k decision trees have been added to the GBDT model; Used as an estimate of the residual; Using a k and β k , jointly adapt to the training data, expressed as: （17）； In the formula, is the regression tree function; S3.3, using the least squares method to fit, a k It is expressed as: （18）； Where a is a sample point in A; β k The optimal value of is determined by: （19）； S3.

4. In each iteration from k=1 to K, the regression tree function used as a leaf node The input space is divided into L non-overlapping regions, represented as , represents the region corresponding to the lth leaf node of the kth tree, and the input space is the optimal input variable set X; In each region, a single constant is predicted as follows: （20）； Where I is an indicator function. When it is established, I=1, otherwise I is equal to 0; yes Estimated residuals in the region The average value of Therefore, in each sub-region In the regression tree, the lth leaf node is solved in the kth iteration. Since the regression tree obtained from equation (20) can predict the , so it is marked as The solution of equation (19) is simplified to: （21）； In the formula, Represents the predicted value of the lth leaf node; S3.

5. In order to supervise the learning rate of the prediction process, all learners are normalized, and the learning rate adjustment factor z is set to adjust the step size of the gradient descent, 0<z≤1. The currently fitted GBDT model is updated in the following way: （22）； In the formula, Represents the prediction function of the GBDT model for X at the kth iteration; S3.

6. When the iteration ends, the final It is expressed as: （23）。 7. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 6 is characterized in that: In S3, the improvement of SGBRT is that a randomization process is introduced in GBDT, that is, a random subset of the training data set is sampled to fit the base learner, and then the updated SGBRT model of the current iteration is calculated based on the separated random subset of the data set, and the steps include: For a given training sample set M, the input space X is divided into L non-overlapping regions ; definition is the loss function, is a set of randomly permuted integers The obtained sequence set, is the i-th sequence; In each sample set of the disjoint region, a sample subset M is randomly selected without replacement s , denoted as , N s is the sample size and N s <N, , Corresponding to x i and i ; Using M s Train the regression tree to obtain a mapping from X to Y; use the random subset as input to generate a new regression tree model for each region, and update equation (21) as: （24）； In the formula, represents the updated solution of equation (19); The final SGBRT model is determined through iteration and is expressed as: （25）。 8. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 7 is characterized in that: In S4, the GWO algorithm is used to further optimize the regression parameters in the SGBRT model, including the number of trees K, the number of disjoint regions L, and the learning rate adjustment factor z.

9. The multivariable wind power curve modeling method based on SGBRT and GWO according to claim 1, characterized in that: In S4, the optimization process of the GWO algorithm is: In S4.1, GWO algorithm, gray wolves are divided into four categories: I, II, III and IV. Wolf I is the chief leader, wolves II and III are subordinates, and the remaining IV wolves are responsible for the balance of relationships within the population. The position of the wolf pack is iteratively updated according to the fitness value of the wolf pack, and its updated position is used to search for the optimal solution in the problem space. During the hunting process, the behavior of gray wolves encircling prey is expressed as: （26）； In the formula, t represents the current number of iterations; , Respectively represent the position vector of the gray wolf at the tth and t+1th iterations; represents the position vector of the prey at the tth iteration; Q and Z are coefficient vectors, respectively expressed as: （27）； Where b1 and b2 are random constant vectors between 0 and 1; q is a convergence factor that decreases linearly from 2 to 0; S4.2, when the prey position is identified, wolves II and III coordinate all wolves to surround the prey under the guidance of wolf I, which is expressed as: （28）； In the formula, , , Respectively represent the position vectors of wolf I, wolf II and wolf III, expressed as: （29）； In the formula, , , They represent the position vectors of wolf I, wolf II and wolf III at the tth iteration respectively; Q1, Q2, Q3 are the coefficient vectors of wolf I, wolf II and wolf III respectively; , , Respectively represent the distances between wolf I, wolf II, and wolf III and other individuals, expressed as: （30）； Where Z1, Z2, and Z3 are random vectors; S4.

3. Constructing the fitness function of GWO , expressed as: （31）； In the formula, RMSE is the root mean square error; MAE is the mean absolute error; R 2 is the coefficient of determination; , , They are RMSE, MAE and R 2 A weighting factor for ; c is a positive constant used to prevent division by zero errors.

Images