Gaussian Process Regression for numerical wind speed prediction enhancement

Renewable Energy 146 (2020) 2112e2123

Contents lists available at ScienceDirect

Renewable Energy journal homepage: www.elsevier.com/locate/renene

Gaussian Process Regression for numerical wind speed prediction enhancement Haoshu Cai, Xiaodong Jia*, Jianshe Feng, Wenzhe Li, Yuan-Ming Hsu, Jay Lee 3NSF I/UCR Center for Intelligent Maintenance Systems, Department of Mechanical Engineering, University of Cincinnati, PO Box 210072, Cincinnati, OH, 45221-0072, USA

a r t i c l e i n f o

a b s t r a c t

Article history: Received 6 January 2019 Received in revised form 27 May 2019 Accepted 4 August 2019 Available online 6 August 2019

This paper studies the application of Multi-Task Gaussian Process (MTGP) regression model to enhance the numerical predictions of wind speed. In the proposed method, a Support Vector Regressor (SVR) is first utilized to fuse the predictions from Numerical Weather Predictors (NWP). The purpose of this regressor is to map the numerical predictions at coarse geographical nodes to the desired turbine location. In subsequent analysis, this SVR prediction output is further enhanced by the MTGP regression model. Based on the validation results on the real-world data from large-scale off-shore wind farm, the prediction accuracies of the NWP are significantly improved at both the short-term horizons (1e6 h ahead) and the long-term horizons (7e24 h ahead) by employing the proposed method. More importantly, the short-term prediction accuracy after enhancement is found comparable or even better than the cutting-edge statistical models for short-term extrapolations. © 2019 Elsevier Ltd. All rights reserved.

Keywords: Wind speed prediction Multi-task Gaussian process Gaussian process regression Support vector machine Time series prediction Forecasting

1. Introduction Driven by the demand of renewable energy, large numbers of wind power generators are erected over the past decade [1e3]. However, one limitation that impedes the further development of wind farm is its high Operation and Maintenance (O&M) cost, which is essentially caused by the uncertainty of wind power production. To better adapt to the inconsistent wind conditions and reduce the costs, wind speed (WS) prediction is identified as one of the key inputs for the wind farm power dispatching and the maintenance planning. In current practices, both the short-term WS prediction for the future 0e6 h and the long-term prediction for future 7e24 h are important inputs to meet the requirements of power grid dispatching, and the power output for each wind turbine are optimized based on the predicted wind conditions [4,5]. For maintenance planning, the short-term WS prediction within future 6 h also serves as a critical input to schedule the maintenance activities for the next day and to minimize the production loss [6e9]. To this end, investigations on the advanced analytics for WS prediction hold great economic value and academic value.

* Corresponding author. 560 Baldwin Hall, University of Cincinnati, PO Box 210072, Cincinnati, OH, 45221, USA. E-mail address: [email protected] (X. Jia). https://doi.org/10.1016/j.renene.2019.08.018 0960-1481/© 2019 Elsevier Ltd. All rights reserved.

WS prediction is intrinsically challenging due to the intermittent fluctuations in WS under intricate meteorological conditions. To better predict the WS at different time horizons, the short-term prediction and long-term prediction are done by different types of data driven models. The short-term prediction is largely based on statistical approaches, which extrapolates the WS series by modeling its time evolution using statistical models. Related examples involve AutoRegressive (AR) and Auto-Regressive Moving Average (ARMA), Artificial Neural Network (ANN) and Kalman Filter (KF) or Unscented Kalman Filter (UKF) that are discussed in Refs. [10e13]. For example, in Ref. [14], linear and non-linear ARMA models are built for 10-min ahead WS prediction. In Ref. [15], an ANN methodology is proposed for super short-term (under 30 s) wind speed prediction, aiming to improve computational efficiency for real-time turbine control. To enhance the performance and the robustness of these prediction algorithms, various combined models are also proposed in the literature. Monfared et al. [12] utilize fuzzy logic and ANN to model WS estimation based on the statistic properties of the input time series. Kani et al. [16] use ANN and Markov Chain to capture patterns in WS time series data. Chen et al. [17] integrate Support Vector Regression (SVR) with KF to realize dynamic state estimation. Santamaría-Bonfil et al. [18,19] utilize SVR model by tuning the model parameters using heuristic algorithms such as

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

Genetic Algorithm (GA) and Particle Swam Optimization (PSO). For the long-term WS prediction, the prediction outputs from the Numerical Weather Predictors (NWP) are generally preferred since the accuracy of the statistical models deteriorates very fast when larger prediction horizon is considered. The prediction results from NWP are usually given at coarse geographical grids and these results may indicate systematic prediction bias on complex terrain. To address these issues, the NWP results are used as reference predictions and regression models are employed to postprocess the NWP outputs. For example, KF is explored as a postprocessing method in Ref. [20] to correct the prediction results from NWP model to avoid systematic bias. Similarly, a practical methodology based on KF is utilized to improve the prediction of NWP in Ref. [21] and the results are validated on two years’ data. ANN is also a useful tool to post-process NWP outputs. Men et al. [22] propose an ensemble of mixture density neural networks to combine the NWP results. Recently, a sequence transfer correction algorithm is also proposed to correct the prediction results from NWP model [23]. In Ref. [24], a Fuzzy System is proposed for NWPbased prediction, and Cuckoo Search (CS) algorithm is utilized to correct the prediction results based on physical laws. In most researches on NWP, prior knowledge about physical and meteorological laws is vital for WS prediction. For example, KF is applied to select forecasting factors such as temperature and air pressure to enhance NWP prediction [25]. In Ref. [26], k-means clustering is combined with ANN to classify NWP data based on meteorological factors. Other techniques such as Support Vector Machine (SVM) and AdaBoost are also combined for NWP clustering and pattern recognition [27]. By summarizing these related researches, it is found that the statistical models for time series extrapolation can give rather satisfactory accuracy in the short-term horizons (1e6 h ahead). However, its accuracy deteriorates very fast when the prediction horizon exceeds 6 h. For long-term prediction (
7 h), the predictions from NWP are necessary to guarantee the prediction accuracy. However, since the NWP results are given at coarse geographical grids and due to the complexity of wind dynamics itself, the direct output from NWP normally has bias comparing with prediction target. This bias in prediction is found more prominent in complex terrain and dynamic wind environment. Although the KF structure in Refs. [20,21] mitigates the prediction bias of NWP to some extent, this method still fails to yield satisfactory results in short-term prediction horizons. Moreover, the KF needs too many prior assumptions, such as the regression coefficients in KF need to be known beforehand, the noise term in KF needs to follow Gaussian distribution and the dynamic process modeled by KF must be linear. To address these challenges and limitations, this work proposes to use the Multiple Task Gaussian Process (MTGP) to post-process the numerical weather predictions. In this work, a novel methodology for NWP enhancement is proposed based on the MTGP model. One major contribution of the proposed method is that it not only enhances the prediction accuracy in the long-term horizons (7e24 h ahead) but also significantly improved the accuracy in the short-term horizons (1e6 h ahead). Especially, the short-term prediction accuracy of the proposed method is found even better than the statistical models that are specifically proposed for short term extrapolations. In the proposed method, spatial correlation of wind speed between the prediction grids of NWP and the turbine nacelle position are first modeled by SVR. Subsequently, the prediction outcome of SVR are further enhanced by MTGP. The superiority of the proposed method is demonstrated by benchmarking with cutting-edge prediction techniques for short term predictions in Ref. [17] and the recent techniques for NWP enhancement in Refs. [20,21]. To the author's knowledge, this is the first time that

2113

MTGP is applied to address wind speed prediction issues. And there are still limited prediction methods in literature that can achieve good accuracy in both short-term and long-term prediction horizons. The rest of this paper is organized as follows. Section 2 illustrates the technical backgrounds. In Section 3, the proposed methodology is presented and described. Section 4 shows the experiment results, the comparison with several benchmarks and the discussions. Finally, conclusions and future work are presented in Section 5.

2. Technical backgrounds 2.1. Standard Gaussian Process Regression Gaussian Process Regression (GPR) is a non-parametric method that can model arbitrary complex system. In most prediction problems, GPR is preferred due to its flexibility to provide the uncertainty representations [28]. GPR models a time series using Gaussian prior that is parameterized by the mean function (MF) mðxÞ and a covariance function (CovF) kðx; x0 Þ as described below:

y ¼ f ðxÞ  NðmðxÞ; kðx; x0 ÞÞ

(1)

In Eq. (1), x and y denote the input and output in the training dataset and f ðxÞ is known as latent variable in the GPR model. In most applications, the mean function mðxÞ in Eq. (1) is set to 0, and CovF kðx; x0 Þ, which describes the similarity between input data points, is the key ingredient in GPR since data points with similar input x are likely to have similar target value y [29]. In the current literature, one of the most frequently-used kernel function squared exponential (SE) is shown below:

kSE ðdÞ ¼ q21 exp



d2

!

2

2q2

(2)

Where d in Eq. (2) are the Euclidean distance between two indexes 0 d ¼ x  x2 . Kv denotes the modified Bessel function. Parameters q1 ,q2 in Eq. (2) are the hyper parameters that need to be optimized. During the model training, the negative log marginalized likelihood (NLML) in Eq. (3) is minimized, so that the hyper-parameters in the kernel matrix K can be estimated.

  1  1 1   NLML ¼  logðpðyjx; qÞÞ ¼  logK þ s2n I  yu K þ s2n I y 2 2 n  logð2pÞ 2 (3) The unknown hyper-parameter q in Eq. (3) is determined by minimizing the NLML. The optimization problem for parameter estimation is written as:

b q ¼ argmin  logðpðyjx; qÞÞ q

(4)

Since the NLML is a convex function, it can be optimized by offthe-shelf optimization algorithms, such as gradient descent. After the model training, the predictive distribution of GPR at testing data point x* can be described as:

   f * x; y; x*  N f * ; covðf * Þ

(5)

2114

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

 1 f * ¼ mðx* Þ þ Kðx; x* Þ Kðx; xÞ þ s2n I ðy  mðxÞÞ  1   Kðx; x* Þ cov f * ¼ Kðx* ; x* Þ  Kðx* ; xÞ Kðx; xÞ þ s2n I

(6)

(7)

Where the f * is the prediction results and covðf * Þ demonstrates the prediction uncertainty. The mean of GPR predictive distribution in Eq. (6) is a linear combination of target variable y in the training set, when the mean function mðxÞ ¼ 0. Under this condition, the mean of predictive distribution can be re-written as:

 1 f * ¼ Kðx; x* Þ Kðx; xÞ þ s2n I y ¼ WGPR y

2.2. Multi-Task Gaussian Process regression MTGP is an extension of the GPR model, and it is described as a special case of standard GPR [31], to deal with the situation when GPR model has multiple outputs. MTGP was originally proposed in Ref. [32], and the superiority of MTGP in the multivariate psychological time-series analysis was demonstrated in Ref. [33]. Another more recent study about using MTGP for battery capacity prediction also presents improved results [30]. In the setting of WS prediction, the input of MTGP is the time indices of the WS series, and the output of MTGP is multiple WS series including the historical WS at turbine nacelle and the reference series from NWP model. The key in MTGP is to recognize the correlation across multiple outputs by using the novel covariance kernel function below:

0

t

t

output series. To ensure that KMTGP is positive semi-definitive, Kc is constructed based on the Cholesky decomposition as below:

2

Step 1: Given training data x and y, the hyper-parameters in the GPR model is obtained by minimizing the NLML in Eq. (3); Step 2: Given the testing time index x , the predictive distribution of GPR is obtained by using Eq. (5) and (6). The hyperparameters in Ref. [Eq. 5~Eq. 7] are obtained from Step 1; Step 3: The mean of the predictive distribution in Eq. (6) is employed as the predicted value, the confidence interval is derived by using the covariance function in Eq. (7). In this study, the error bound is as f * ±2covðf * Þ.

(9)

where l; l 2f1; 2; …; mg represent the indices of series and there are m series in total, kc and kt in Eq. (9) model the correlation across the multiple outputs and the covariance for one series respectively. x; x0 denote the time indices for task l and l0 . Based on Eq. (9), the kernel matrix of MTGP can be constructed as:

(10)

where 5 is the Kronecker product. qc and qt are the hyperparameters in the kernel matrix. In Eq. (10), Kc is a m  m similarity matrix that models the correlation or similarity across multiple series, Kt is a nt  nt symmetric matrix that models the covariance across the time indices for the t-th series , nt represents the number of time indices for the t-th serie. Therefore, KMTGP is a P t  P t  n  n matrix that captures the similarity across multiple

(8)

Where WGPR is the weighting matrix for the standard GPR. The mean function of GPR mðxÞ is normally set to be 0 for trendfree time series. GPR is known as non-parametric approach which can be employed to model time series or systems with arbitrary complexity when provided with sufficient data. A non-zero mean function is normally employed when a clear trend is observed from the time series or there is a sound assumption of the trend term. Like in Ref. [30], an exponential trend term is employed as the mean function to better extrapolate the degradation trajectory of the battery cell in the long term. In the present study, the mean function is set to 0 since we did not see a consistent trend term of wind speed in the prediction horizon. When using GPR to make prediction, there are several general steps to follow:

 0  0 0 0 kMTGP x; x ; l; l ¼ kc l; l  kt ðx; x Þ

KMTGP ðX; L; qc ; qt Þ ¼ Kc ðL; qc Þ5Kt ðX; qt Þ

Kc ¼ LLu

6 qcð1;1Þ 6 qcð2;1Þ ¼6 6 « 4

qcð1;2Þ qcð2;2Þ «

qcðm;1Þ qcðm;2Þ

3

/ … 1 …

qcð1;mÞ 7 qcð2;mÞ 7 7 «

qcðm;mÞ

7 5

(11)

where L is a lower triangular matrix. The elements in Kc represent the similarity level between each pair of the reference series. According to the description in Refs. [32,33], these elements in Kc can be interpreted as correlation coefficients.

3. Methodology 3.1. Using MTGP for NWP enhancement When using MTGP for time series prediction, all the training and testing procedures are the same with traditional GPR except the construction of the kernel matrix. To better illustrate the kernel matrix in MTGP, we use the case of two output tasks as an example. In this scenario, the MTGP model can be constructed as:

2

3 0 2 y r 4 5 ¼ N @0; 4 qrr Kðx  r ; xr Þ yp qpr K xp ; xr ¼ Nð0; KMTGP ðx; xÞÞ





31

qrp Kxr ; xp  5A qpp K xp ; xp (12)

where yr 2Rr1 and yp 2R p1 are two outputs with different dimensionality, xr and xp are the models input for the two output tasks, the dimensionality of xr and xp should be the same. In 1D case or time series prediction, xr and xp are the time indices of the two time series.qrr , qrp , qpr and qpp denote the correlation coefficients of two output series, these four coefficient are treated as part of the hyper-parameters in MTGP and they are obtained by optimizing the NLML in the model training. It is also important to note that qrp ¼ qpr is always valid due to Eq. (11). Application of MTGP for NWP enhancement is illustrated in Fig. 1. In Fig. 1, yr serves as the reference wind speed series that is given by the NWP and its time indices of yr is written as xr ¼ fTtk ; …; Tt ; …; Ttþ24 g. yp in Eq. (12) corresponds to the measured wind speed series at the turbine nacelle which is going to be extrapolated to the future 24 h, the time indices of yp is as xp ¼ fTtk ; …; Tt g, where k is the length of time block for model construction. When applying KMTGP to the model, it is recommended to min-max normalized the time indices xp into range ½0; 1, and xr is normalized based on the same minimum and maximum of xp . It is just a numerical treatment to avoid very small model parameters and potential under flow. By using yr , yp , xr and xp as training data, a MTGP model is constructed and a better prediction at time x*p ¼ fTtþ1 ; …; Ttþ24 g can be obtained. The predictive distribution can be simply obtained following Eq. (5)e (7). Like the standard GPR, the prediction results of MTGP can be

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

2115

Fig. 1. MTGP enhancement for reference-based prediction.

also interpreted as a linear combination of the historical observations, which can be written as:

* fp

¼



KMTGP x; x*p

Where

2 3 2 3  1 y y r 5 2 4 KMTGP ðx; xÞ þ sn I ¼ WMTGP 4 r 5 yp yp

2





3

q Kðx ; x Þ qrp Kxr ; xp  5 KMTGP ðx; xÞ ¼ 4 rr  r r  qpr K xp ; xr qpp K xp ; xp      KMTGP x; x*p ¼ qpr K x*p ; xr



qpp K x*p ; xp

ðrþpÞðrþpÞ



This is achieved by obtaining the optimal qpp and qpr in the model training phase. Therefore, the model is expected to be superior than NWP in long-term horizons and to be superior than statistical models in short-term horizons. As a summary of the discussion, the algorithm proposed for NWP enhancement based on MTGP is described as below. It is important to note that the proposed method requires to re-train the MTGP model at each prediction step. This implies that the proposed model has no bias known as seasonal effects, since a new model is derived at each prediction step by utilizing the data in the short past only. Algorithm 1. Using MTGP for NWP Output Enhancement

24ðrþpÞ

(13) From Eq. (13), one can easily find that the prediction output of

At certain time step , Step 1:

Initialize

,

,

,

,

and the time block length

, is the NWP prediction at time nacelle during time

,

as:

, , is the measured wind speed at turbine

Step 2:

Construct the MTGP as Eq. Error! Reference source not found. and obtain the optimized hyper-parameters , , and other hyper-parameters in the kernel function. The optimal hyper-parameters are obtained by minimizing the NLML.

Step 3:

Obtain the predictive distribution of MTGP. The prediction mean is described as Eq. Error! Reference source not found.).

Step 4:

Propagate to time

and repeat Step 1~Step 3.

the MTGP model is simply a linear combination of yr and yp . Consequently, MTGP is employed as a novel approach to further enhance the prediction accuracy of NWP in this investigation. Due to the fact that the prediction result of GPR and MTGP is given as a linear combination of yr and yp , the MTGP model enhances the NWP results both in the short-term and long-term accuracy. The presence of qpp Kðxp ; xp Þ term in the kernel matrix mainly contributes to the short-term prediction accuracy, since the qpp is more dominant than qpr in the short-term horizons. The prediction result considers yr as a reference in the long-term extrapolations, which is why the prediction accuracy does not deteriorate significantly in the long-term horizons. More importantly, the algorithms can achieve enhanced prediction accuracy in both short-term and long-term horizons mainly because MTGP automatically decides the optimal trade-off between the NWP output and the extrapolation of measured wind speed at turbine nacelle.

3.2. The proposed methodology and implementation The goal of this research is to predict WS in future 24 h. The given data in this research involves the WS at the turbine nacelle collected by the Supervisory Control and Data Acquisition (SCADA) system and the weather forecast data from the NWP model. The NWP data is given as the average WS within 1 h at each forecast gird from three different heights, 10m, 50m and 90m. In this research, the NWP data from today and one day before is used to build the model. Therefore, at any specific point of time, the available data includes the SCADA data up till now, and the NWP data updated today and 1 day ago. Spatially, the NWP data includes the WS prediction at 9 grid nodes in Fig. 2, which has rather coarse spatial resolution of roughly 9.5 km. The Wind Turbine (WT) position is located within the area covered by these grid nodes, as shown in Fig. 2, and the heights of the turbines remain unknown.

2116

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

RMSEh ¼

vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi uP  2 0 u m t i¼1 y*p;h  yp;h m

  0  m y*  y 1 X  p;h p;h  MAPEh ¼   0  m i¼1  yp;h

Fig. 2. Wt position and weather forecast grid position.

The proposed methodology is illustrated in Fig. 3. In Fig. 3, the time resolution of NWP data is 1 h. As mentioned before, the NWP data reported from today and one day before is utilized. Therefore, the dimensionality of numerical weather predictions is 3  9  2 ¼ 54. As shown in Fig. 3, a time block with length k is needs to establish the SVR model between all 54 numerical wind predictions and the measured wind speed series. The purpose of this SVR model is to give a reference prediction series yr in the future 24 h. This prediction is merely based on the NWP outputs and will be subsequently enhanced by the MTGP model. In the MTGP enhancement step, the detailed procedures for model construction and making prediction are described in Algorithm 1. It is important to highlight that the SVR model and the MTGP model are re-trained at each time step by using the historical data for model construction. Therefore, the seasonal effect of the wind speed distribution is not a concern in the present model, because the prediction is made based on the predication output of NWP and the extrapolation of wind series in the recent past. To calibrate the performance of the propose method, Root Mean Square Error (RMSE) and Mean Absolute Percentage Error (MAPE) are used as criteria to evaluate the prediction accuracy. Suppose the 0 WS data from current time point to the next 24 h is given by yp ¼ 0

0

0

0

fyp;1 ; yp;2 ; :::yp;h ; :::yp;24 g, the predicted WS data is :::y*p;h ;:::y*p;24 g,

y*p

¼

fy*p;1 ; y*p;2 ;

RMSE and MAPE at a certain prediction horizon h are

calculated as follows:

(14)

(15)

Where m represents the number of prediction steps, h ¼ 1; …; 24 denotes the prediction horizon. To better demonstrate the improvements made by the proposed methodology, the methods listed in Table 1 will be benchmarked. In Table 1, Best NWP model uses RMSE to select one series of NWP data with the highest forecast accuracy. SVR-NWP model uses an SVR model to fuse the NWP data for prediction as Step 1 in Fig. 3. Similarly, GPR-NWP model uses a GPR model to fuse the NWP data for prediction. SVR þ UKF focuses on short-term prediction and it is a dynamic methodology that is proposed in Ref. [17]. GPR model extrapolates historical WS data at each time step, which focuses on short-term prediction as well. The employed GPR model is implemented with SE kernel function. To compare with other peer algorithms for NWP enhancement, KF structure that is discussed in Refs. [20,21] are implemented and benchmarked. In their discussions, the bias of the NWP is modeled as a high-order polynomial: *3 et ¼ c0;t þ c1;t ,x*t þ c2;t ,x*2 t þ c3;t ,xt þ vt

(16)

where et is a scalar bias of NWP at time t, c0;t , c1;t , c2;t and c3;t are the polynomial coefficients. vt is Gaussian process noise. The above equation is implemented in a KF structure as shown below:

i h *3 ct ¼ ct1 þ wt ; et ¼ 1; x*t ; x*2 t ; xt ,ct þ vt

(17)

where ct ¼ ½c0;t ; c1;t ; c2;t ; c3;t u . To summarize, a list of benchmarking algorithms is tabulated in Table 1. And the prediction accuracies the algorithms in Table 1 are benchmarked based on the real-word data in the next section. 4. Results and discussions The performance of the proposed method is validated based on the off-shore wind farm data collected within half a year. This dataset under study includes the WS series at turbine nacelle collected by the SCADA and the NWP forecast data that is described above. For KF ensembled methods in Ref. [20], the data from January and February 2017 is used for cross-validation and model training, while data from April and September 2017 is used for

Fig. 3. Statement of the WS prediction problem.

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

2117

Table 1 List of benchmarking algorithms. Short-term prediction (1e6 h ahead)

Long-term prediction (6e24 h ahead)

Short and long-term prediction (1e24 h ahead)

SVR þ UKF [17] GPR

Best NWP þ KF [21] SVR-NWP þ KF [21] SVR-NWP Best NWP GPR-NWP

SVR-NWP þ MTGP (The proposed method)

Table 2 The parameters and distributions of WS. Turbine #1

Turbine #2

April

September

April

September

Max

Min

Mean

STD

Max

Min

Mean

STD

Max

Min

Mean

STD

Max

Min

Mean

STD

26.73

0.67

7.55

4.31

19.65

0.57

6.42

3.34

25.62

0.60

7.12

3.81

20.53

0.76

6.37

3.33

Distribution of WS

80

60

60

60

40 20 0

Count

80

40 20

5

0

10 15 20 25

WS(m/s)

Distribution of WS

60

Count

Distribution of WS

80

Count

Count

Distribution of WS

40

20

20 5

10

15

0

20

40

5

WS(m/s)

0

10 15 20 25

5

WS(m/s)

10

15

20

WS(m/s)

Table 3 The information of NWP on grid #1, April, reported from today. Height 10m

Height 50m

Max

Min

Mean

STD

17.30

0.03

6.22

3.40

Distribution of NWP

Height 90m

Max

Min

Mean

STD

19.78

0.25

7.22

3.65

Distribution of NWP

Max

Min

Mean

STD

20.89

0.12

7.70

3.82

Distribution of NWP

20 0

60

Count

40

Count

Count

60 40 20

5

10

15

WS(m/s)

20

0

40 20

5

10

15

WS(m/s)

model performance testing. For other methods which don't need pre-training, the data from April and September is used for model testing and result benchmarking. At the beginning of the analysis, all the time series from NWP and anemometer measurements at turbine nacelle are synchronized and pre-processed to have a 1-h interval. To demonstrate the advantages of the proposed methodology, the prediction model is validated on two different turbines. Table 2 presents the information of WS recorded on the two target turbines. Table 3 and Table 4 present the information of NWP by taking two examples from Grid #1 and Grid 9 on April and September respectively. The prediction results are benchmarked in Table 5 ~ Table 8. The SVR model employed in benchmarks and the proposed method is built with a gaussian kernel, and the hyper-

20

0

5

10

15

20

25

WS(m/s)

parameters are validated by gird-searching. Tables 5 and 6 show the prediction results of WT #1 in April and September. Table 7 and Table 8 show the prediction results of WT #2 during the same months. The best prediction accuracies at each prediction horizon are highlighted in bold character. Generally, the proposed model SVR-NWP þ MTGP yields the best results comparing with all others. Comparing the prediction accuracies of different methods at different prediction horizons, the following findings are highlighted: 1) GPR model and SVR þ UKF model proposed in Ref. [17] demonstrate improved prediction accuracy in short-term

2118

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

Table 4 The information of NWP on grid #9, September, reported from yesterday. Height 10m

Height 50m

Max

Min

Mean

STD

14.73

0.18

6.18

2.67

Distribution of NWP

Height 90m

Max

Min

Mean

STD

18.63

0.41

7.31

3.19

Distribution of NWP

Max

Min

Mean

STD

21.43

0.42

8.03

3.67

Distribution of NWP

80

60

60

40 20 0

Count

Count

Count

60 40

20

20

5

10

0

15

40

0 5

WS(m/s)

10

15

20

5

10

15

20

WS(m/s)

WS(m/s)

Table 5 Turbine #1 comparison result of April. Methods

Predict Hours 1- hour ahead

SVR þ UKF [17] GPR Best NWP SVR-NWP SVR-NWP þ KF [21] Best NWP þ KF [21] GPR-NWP SVR-NWP þ MTGP (The proposed method)

4-h ahead

6-h ahead

12-h ahead

18-h ahead

24-h ahead

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

21.578 17.404 50.810 32.462 30.179 37.949 30.599 18.891

1.320 1.175 3.729 2.404 2.069 2.689 1.776 1.281

42.362 46.205 57.339 45.674 43.264 49.353 42.188 41.312

2.960 2.988 3.951 3.302 2.747 3.450 2.528 2.481

46.781 53.678 59.534 49.068 46.452 50.793 45.492 44.645

3.504 3.853 4.009 3.571 2.896 3.364 2.755 2.669

56.992 65.720 60.912 52.059 50.164 56.417 47.614 49.831

4.516 5.102 3.993 3.994 3.109 3.834 3.043 2.872

57.390 68.381 57.215 52.630 51.725 58.062 49.188 52.026

4.780 5.437 3.874 4.225 3.274 3.915 3.237 2.991

59.211 69.144 48.689 55.041 55.526 58.855 52.760 51.488

5.192 5.673 3.551 4.404 3.346 3.821 3.325 3.013

Table 6 Turbine #1 comparison result of September. Methods

SVR þ UKF [17] GPR Best NWP SVR-NWP SVR-NWP þ KF [21] Best NWP þ KF [21] GPR-NWP SVR-NWP þ MTGP (The proposed method)

Predict Hours 1- hour ahead

4-h ahead

12-h ahead

18-h ahead

24-h ahead

MAPE

RMSE

MAPE

RMSE

MAPE

6-h ahead RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

24.493 19.640 61.311 34.850 35.316 37.598 36.045 22.588

1.437 1.313 3.448 2.267 2.585 2.670 2.217 1.442

40.850 40.056 64.257 43.685 46.288 44.393 46.915 42.071

2.656 2.405 3.512 2.709 2.826 3.147 2.788 2.577

46.770 47.404 63.280 45.656 49.599 46.075 48.234 45.660

3.342 2.826 3.488 2.807 2.940 3.321 2.875 2.780

57.907 58.445 60.932 48.943 56.513 48.928 51.496 48.101

4.440 3.476 3.474 2.959 3.206 3.365 3.007 2.853

64.258 67.045 64.637 49.246 58.604 48.732 51.627 48.635

5.051 3.939 3.561 2.999 3.312 3.386 3.042 2.919

69.382 71.048 59.814 50.100 58.544 53.402 53.922 50.383

5.488 4.269 3.332 3.077 3.387 3.593 3.132 3.109

Table 7 Turbine #2 comparison result of April. Methods

Predict Hours 1- hour ahead

SVR þ UKF [17] GPR Best NWP SVR-NWP SVR-NWP þ KF [21] Best NWP þ KF [21] GPR-NWP SVR-NWP þ MTGP (The proposed method)

4-h ahead

6-h ahead

12-h ahead

18-h ahead

24-h ahead

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

20.522 16.947 54.510 32.938 31.651 36.036 30.490 18.317

1.284 1.169 3.511 2.256 1.900 2.563 1.799 1.240

42.853 44.280 59.528 43.780 43.678 49.136 44.029 38.614

2.597 2.830 3.663 3.067 2.577 3.424 2.574 2.274

48.267 53.271 61.808 46.679 46.279 50.643 46.388 41.994

3.030 3.598 3.727 3.328 2.750 3.379 2.776 2.466

58.061 64.750 64.229 50.568 49.773 54.012 48.722 46.547

3.993 4.430 3.722 3.682 2.962 3.545 2.986 2.656

62.880 66.266 62.155 50.279 51.621 54.849 50.021 47.635

4.487 4.813 3.659 3.824 3.057 3.508 3.100 2.755

60.653 64.459 52.734 52.941 53.239 58.816 52.071 48.145

4.612 4.969 3.380 3.978 3.115 3.899 3.147 2.797

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

2119

Table 8 Turbine #2 comparison result of September. Methods

Predict Hours

SVR þ UKF [17] GPR Best NWP SVR-NWP SVR-NWP þ KF [21] Best NWP þ KF [21] GPR-NWP SVR-NWP þ MTGP (The proposed method)

1- hour ahead

4-h ahead

6-h ahead

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

MAPE

RMSE

23.349 17.517 59.158 33.385 33.409 34.456 33.112 21.661

1.387 1.250 3.441 2.343 2.268 2.539 2.164 1.452

41.875 35.750 60.907 41.175 42.550 41.002 43.583 39.950

2.618 2.343 3.460 2.757 2.787 3.071 2.793 2.564

45.514 41.984 62.713 43.057 45.107 43.953 45.490 44.061

3.124 2.732 3.609 2.855 2.921 3.334 2.912 2.848

55.323 54.333 60.261 45.717 50.296 46.819 50.166 47.951

3.994 3.751 3.532 2.954 3.139 3.460 3.083 3.062

65.124 63.085 64.388 46.806 50.496 46.340 51.282 47.378

4.764 4.437 3.649 3.018 3.210 3.465 3.166 3.041

71.706 73.565 59.738 48.089 50.851 51.382 52.501 48.272

5.212 5.071 3.461 3.098 3.271 3.614 3.236 3.146

horizons (1e6 h ahead). However, their extrapolation accuracies deteriorate very fast at long-term horizons (7e24 h ahead); 2) GPR-NWP, Best NWP and SVR-NWP demonstrate good accuracy at long-term horizons. However, their accuracies in short-term are not comparable to the statistical models; 3) Best NWP þ KF and SVR-NWP þ KF demonstrate enhanced prediction accuracy comparing with Best NWP and SVR-NWP in short-term horizons, especially 1e3 h ahead. This finding indicates that the KF structure in Refs. [20,21] can effectively

6

(a)

short term 60

short term

50

4

MAPE(%)

RMSE(m/s)

4.5

3.5 3

2 1.5 0

5

10

15

20

Turbine #1, April 40 long term 30

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

long term

2.5

1

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

20

10

25

0

5

5.5

15

20

25

80

(d) Turbine #1, September

5

10

Predicted Hour(h)

Predicted Hour(h)

(c)

24-h ahead

70

Turbine #1, April

5

18-h ahead

enhance the prediction accuracy as expected. However, the performance of such post-processing steps in long-term prediction horizon is quite unstable. In some occasions, it makes the prediction accuracy even worse; 4) The proposed method, SVR-NWP þ MTGP, demonstrates excellent accuracy in both short-term and long-term horizons. Its prediction accuracies in short-term horizons are found comparable to SVR þ UKF model and even better in some occasions. More importantly, the long-term prediction accuracy of

(b)

5.5

12-h ahead

short term

Turbine #1, September

70

short term

4.5

60

MAPE(%)

RMSE(m/s)

4 3.5 3 2.5

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

long term

2 1.5 1

0

5

10

15

Predicted Hour(h)

20

50 long term

40

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

30 20

25

10

0

5

10

15

20

Predicted Hour(h)

Fig. 4. Benchmarking of Turbine #1 prediction. (a) RMSE in April; (b) MAPE in April; (c) RMSE in September; (d) MAPE in September.

25

2120

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

the proposed method is better than Best NWP, SVR-NWP, Best NWP þ KF and SVR-NWP þ KF. In addition, the improvements made by the proposed method is consistent over different turbine locations and different month. The superiority of the proposed method is better explained in Fig. 4 and Fig. 5 by comparing the RMSE and MAPE of different methods. The results in Figs. 4 and 5 demonstrate the validation results on two different wind turbines at two different months. In both short-term and long-term prediction horizons, the proposed method gives the best accuracy in term of both RMSE and MAPE. It is highlighted that the short-term prediction performance of the proposed method is comparable with the state-of-art approach SVR þ UKF in Ref. [17], the long-term prediction accuracy of the proposed method is superior than the KF structure that is presented in Ref. [21] for NWP post-processing. Fig. 6 shows detailed results of SVR-NWP þ MTGP model, which is randomly selected from the testing dataset. Several findings are highlighted below: (1) Generally, the proposed model fits the ground truth well. (2) Short-term horizon leads to smaller error bounds. For long-term horizon, the proposed model is capable to describe the overall trends of WS data. (3) In most cases, the true value of WS data falls into the predicted error bounds, except the

5

(a)

wind gust around time point 300e350. However, for 1-h ahead prediction, the proposed model is capable to predict the wind gust very well. Fig. 7 compares the error distributions of the benchmarking models and the proposed model in two months. Several key points for discussion in Fig. 7 are displayed as follows: (1) The SVR þ UKF and the proposed method give the smallest error in 1-h ahead prediction. And the proposed method is slightly better than SVR þ UKF; (2) the proposed method also gives the best prediction at 6-h ahead prediction. The prediction error of the SVR þ UKF method starts to increase at this prediction horizon; (3) At 12-h ahead prediction and 24-h ahead prediction, the SVR-NWP and SVR-NWP þ MTGP give the smallest error. (4) Comparing SVR þ UKF, Best NWP and the proposed method, one can easily find that the proposed method keeps the advantage of SVR þ UKF in the short-term horizons and the advantage of NWP in the long-term horizons; (5) Comparing Best NWP, SVR þ NWP and SVRNWP þ MTGP, one can find that the MTGP mainly reduces the prediction error in the short-term horizons. At long-term horizons, the prediction results are slightly better than SVR-NWP and significantly better than NWP; (6) By comparing with the recent models in literature SVR þ UKF [17] and NWP þ KF [21], the proposed method gives the best overall prediction accuracy in the 24 h

70

(b)

short term

Turbine #2, April

4.5

60

short term

50

3.5

MAPE(%)

RMSE(m/s)

4

3 2.5

1.5 1

0

5

10

15

20

long term 30

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

long term 2

Turbine #2, April 40

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

20

10

25

0

5

5.5

(c)

10

15

20

25

Predicted Hour(h)

Predicted Hour(h) 80

(d)

Turbine #2, September

short term

5

Turbine #2, September

70

short term

4.5

60

MAPE(%)

RMSE(m/s)

4 3.5 3 2.5

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

long term

2 1.5 1

0

5

10

15

Predicted Hour(h)

20

50 long term

40

SVR+UKF[17] GPR Best NWP SVR-NWP SVR-NWP+KF[21] Best NWP+KF[21] GPR-NWP SVR-NWP+MTGP

30 20

25

10

0

5

10

15

20

Predicted Hour(h)

Fig. 5. Benchmarking of Turbine #2 prediction. (a) RMSE in April; (b) MAPE in April; (c) RMSE in September; (d) MAPE in September.

25

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

2121

Error Bound Predicted Value True Value

1 hour Ahead

WS(m/s)

15 10 5 0 50

100

150

200

(a)

250

300

350

400

450

Time(Hours) 6 hours Ahead

20

Error Bound Predicted Value True Value

WS(m/s)

15 10 5 0 50

100

150

200

250

(b)

300

350

400

450

Time(Hours) 20

Error Bound Predicted Value True Value

12 hours Ahead

WS(m/s)

15 10 5 0 50

100

150

200

250

(c)

300

350

400

450

Time(Hours)

WS(m/s)

Error Bound Predicted Value True Value

24 hours Ahead

15 10 5 0 50

100

150

200

(d)

250

300

350

400

450

Time(Hours)

Fig. 6. Detailed Results and Error Bounds of SVR-NWP þ MTGP model (a) 1-h Ahead Prediction; (b) 6-h Ahead Prediction; (c) 12-h Ahead Prediction; (d) 24-h Ahead Prediction.

ahead prediction. Finally, the execution time of the proposed model is discussed. The proposed methodology and the benchmarking methods are run on a PC with RAM 32 GB, CPU 3.50 GHz, Windows 10 Enterprise. At each time point, WS data and NWP data of the last 120 h are included in all of the models testing procedure. The proposed model is run 10 times on the data from April and September. Fig. 8 demonstrates the average execution time and the accuracy

of different time block lengths. The average execution time refers to the average prediction time at each time point. The time block length is donated to the number of hours during which past WS and NWP data is taken into consideration, as shown in Fig. 3. It indicates that the execution time grows as the length of time block enlarges. Meanwhile, the prediction accuracies RMSE and MAPE are improved as more WS and NWP data is taken into the model. However, when the time block length grows to 96 h, the average

2122

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

1-hour Ahead Count

SVR+ UKF [17]

Count

Best NWP

Count

Best NWP+ KF [21]

SVRNWP+ MTGP

Count

Count

SVRNWP

6-hour Ahead

12-hour Ahead

24-hour Ahead

200

100

100

100

100

50

50

50

0 -15 -10 -5 0 5 10 15

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

0 -20-15-10 -5 0 5 10 15 20

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

200

100

100

100

100

50

50

50

0 -15 -10 -5 0 5 10 15

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

0 -20-15-10 -5 0 5 10 15 20

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

200

100

100

100

100

50

50

50

0 -15 -10 -5 0 5 10 15

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

0 -20-15-10 -5 0 5 10 15 20

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

200

100

100

100

100

50

50

50

0 -15 -10 -5 0 5 10 15

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

0 -20-15-10 -5 0 5 10 15 20

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

200

100

100

100

100

50

50

50

0 -15 -10 -5 0 5 10 15

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

0 -20-15-10 -5 0 5 10 15 20

0 -20-15-10 -5 0 5 10 15 20

Residue

Residue

Fig. 7. Probabilistic histograms of residues with the distribution mean and 2-s interval.

Average Execution Time(s)

6 5 4 3 2 1 0 24

48

(a)

72

96

120

Time Block Length(h) 80

5

70

MAPE(%)

RMSE

4

3 24 hours 48 hours 72 hours 96 hours 120 hours

2

1

(b)

60 50 40 30 20

0

5

10

15

Predicted Hour(h)

20

24 hours 48 hours 72 hours 96 hours 120 hours

25

(c)

0

5

10

15

20

25

Predicted Hour(h)

Fig. 8. Average Execution Time and Prediction Accuracy of Time Block Length. (a) Average Execution Time vs. Time Block Length (b) RMSE of Different Time Block Length (c) MAPE of Different Time Block Length.

H. Cai et al. / Renewable Energy 146 (2020) 2112e2123

execution time still grows while RMSE and MAPE converge to a stable level. Therefore, a time block of 96 h is recommended to achieve the optimal prediction performance and to save execution time. 5. Conclusion In this paper, a novel WS prediction method is proposed. The effectiveness and the superiority of the proposed method are validated on a dataset collected from an off-shore wind farm. The final results suggest following conclusions. (1) The proposed method can be effectively employed to improve the prediction accuracy of the numerical weather prediction; (2) the proposed method carries both the advantages of time series extrapolation method for short-term prediction and the advantages of NWP in long-term prediction horizon. (3) The proposed method reports improved prediction accuracy comparing with the recently proposed models of SVR þ UKF [17] and NWP þ KF [21]. In future works, the proposed method will be further studied and the sparse Gaussian process methods will be explored to further boost the computational efficiency. References [1] I. Colak, S. Sagiroglu, M. Yesilbudak, Data mining and wind power prediction: a literature review, Renew. Energy 46 (2012) 241e247. [2] J. Jung, R.P. Broadwater, Current status and future advances for wind speed and power forecasting, Renew. Sustain. Energy Rev. 31 (2014/03/01/2014) 762e777. [3] X. Jia, C. Jin, M. Buzza, Y. Di, D. Siegel, J. Lee, A deviation based assessment methodology for multiple machine health patterns classification and fault detection, Mech. Syst. Signal Process. 99 (2018) 244e261. [4] J. Jin, D. Zhou, P. Zhou, S. Qian, M. Zhang, Dispatching strategies for coordinating environmental awareness and risk perception in wind power integrated system, Energy 106 (2016) 453e463. [5] X. Zhang, W. Cai, Z. Gan, Optimal dispatching strategies of active power for DFIG wind farm based on GA algorithm, in: Control and Decision Conference (CCDC), 2016 Chinese, 2016, pp. 6094e6099. [6] P. Eecen, H. Braam, L. Rademakers, T. Obdam, Estimating costs of operations and maintenance of offshore wind farms, in: European Wind Energy Conference and Exhibition, Milan, Italy, 2007.
cs, G. Erdo €s, L. Monostori, Z.J. Viharos, Scheduling the maintenance of [7] A. Kova wind farms for minimizing production loss, IFAC Proceedings Volumes 44 (2011) 14802e14807. }s, Z.J. Viharos, L. Monostori, A system for the detailed [8] A. Kovacs, G. Erdo scheduling of wind farm maintenance, CIRP Ann. - Manuf. Technol. 60 (2011) 497e501. [9] X. Jia, C. Jin, M. Buzza, W. Wang, J. Lee, Wind turbine performance degradation assessment based on a novel similarity metric for machine performance curves, Renew. Energy 99 (2016) 1191e1201. [10] H. Liu, H.-q. Tian, Y.-f. Li, Comparison of two new ARIMA-ANN and ARIMAKalman hybrid methods for wind speed prediction, Appl. Energy 98 (2012) 415e424. [11] O.B. Shukur, M.H. Lee, Daily wind speed forecasting through hybrid KF-ANN model based on ARIMA, Renew. Energy 76 (2015) 637e647, 2015/04/01. [12] M. Monfared, H. Rastegar, H.M. Kojabadi, A new strategy for wind speed forecasting using artificial intelligent methods, Renew. Energy 34 (2009)

2123

845e848, 2009/03/01. [13] X. Jia, Y. Di, J. Feng, Q. Yang, H. Dai, J. Lee, Adaptive virtual metrology for semiconductor chemical mechanical planarization process using GMDH-type polynomial neural networks, J. Process Control 62 (2018) 44e54, 2018/02/01. [14] M. Lydia, S. Suresh Kumar, A. Immanuel Selvakumar, G. Edwin Prem Kumar, Linear and non-linear autoregressive models for short-term wind speed forecasting, Energy Convers. Manag. 112 (2016) 115e124, 2016/03/15. [15] J. Luna, S. Gros, J. Geisler, O. Falkenberg, R. Noga, A. Schild, Super-short term wind speed prediction based on artificial neural networks for wind turbine Control applications, in: IECON 2018 - 44th Annual Conference of the IEEE Industrial Electronics Society, 2018, pp. 1952e1957. [16] S.A. Pourmousavi Kani, M.M. Ardehali, Very short-term wind speed prediction: a new artificial neural networkeMarkov chain model, Energy Convers. Manag. 52 (2011) 738e745, 2011/01/01. [17] K. Chen, J. Yu, Short-term wind speed prediction using an unscented Kalman filter based state-space support vector regression approach, Appl. Energy 113 (2014) 690e705. [18] G. Santamaría-Bonfil, A. Reyes-Ballesteros, C. Gershenson, Wind speed forecasting for wind farms: a method based on support vector regression, Renew. Energy 85 (2016) 790e809, 2016/01/01.

rez-Bellido, A. Portilla-Figueras, [19] S. Salcedo-Sanz, E.G. Ortiz-García, A.M. Pe L. Prieto, Short term wind speed prediction based on evolutionary support vector regression algorithms, Expert Syst. Appl. 38 (2011) 4052e4057, 2011/ 04/01. [20] P. Louka, G. Galanis, N. Siebert, G. Kariniotakis, P. Katsafados, I. Pytharoulis, et al., Improvements in wind speed forecasts for wind power prediction purposes using Kalman filtering, J. Wind Eng. Ind. Aerodyn. 96 (2008) 2348e2362. [21] F. Cassola, M. Burlando, Wind speed and wind energy forecast through Kalman filtering of Numerical Weather Prediction model output, Appl. Energy 99 (2012) 154e166, 2012/11/01. [22] Z. Men, E. Yee, F.-S. Lien, D. Wen, Y. Chen, Short-term wind speed and power forecasting using an ensemble of mixture density neural networks, Renew. Energy 87 (2016) 203e211, 2016/03/01. [23] H. Wang, S. Han, Y. Liu, J. Yan, L. Li, Sequence transfer correction algorithm for numerical weather prediction wind speed and its application in a wind power forecasting system, Appl. Energy 237 (2019) 1e10, 2019/03/01. [24] J. Zhao, Z.-H. Guo, Z.-Y. Su, Z.-Y. Zhao, X. Xiao, F. Liu, An improved multi-step forecasting model based on WRF ensembles and creative fuzzy systems for wind speed, Appl. Energy 162 (2016) 808e826, 2016/01/15. [25] S. Hua, S. Wang, S. Jin, S. Feng, B. Wang, Wind speed optimisation method of numerical prediction for wind farm based on Kalman filter method, J. Eng. (2017) 1146e1149, 2017. [26] L. Dong, L. Wang, S.F. Khahro, S. Gao, X. Liao, Wind power day-ahead prediction with cluster analysis of NWP, Renew. Sustain. Energy Rev. 60 (2016) 1206e1212, 2016/07/01. € [27] S. Buhan, Y. Ozkazanç, I. Çadırcı, Wind pattern recognition and reference wind mast data correlations with NWP for improved wind-electric power forecasts, IEEE Transactions on Industrial Informatics 12 (2016) 991e1004. [28] C.E. Rasmussen, Gaussian processes in machine learning, in: O. Bousquet, €tsch (Eds.), Advanced Lectures on Machine Learning: ML U. von Luxburg, G. Ra Summer Schools 2003, Canberra, Australia, February 2 - 14, 2003, Tübingen, Germany, August 4 - 16, 2003, Revised Lectures, Springer Berlin Heidelberg, Berlin, Heidelberg, 2004, pp. 63e71. [29] C.E. Rasmussen, Gaussian Processes for Machine Learning, 2006. [30] R.R. Richardson, M.A. Osborne, D.A. Howey, Gaussian process regression for forecasting battery state of health, J. Power Sources 357 (2017) 209e219. [31] H. Liu, J. Cai, Y.-S. Ong, Remarks on multi-output Gaussian process regression, Knowl. Based Syst. 144 (2018) 102e121. [32] E.V. Bonilla, K.M. Chai, C. Williams, Multi-task Gaussian process prediction, in: Advances in Neural Information Processing Systems, 2008, pp. 153e160. [33] R. Dürichen, M.A. Pimentel, L. Clifton, A. Schweikard, D.A. Clifton, Multitask Gaussian processes for multivariate physiological time-series analysis, IEEE (Inst. Electr. Electron. Eng.) Trans. Biomed. Eng. 62 (2015) 314e322.