Wavelet-Gaussian Process Regression Model for Regression Daily Solar Radiation in Ghardaia, Algeria

2.1 Study region and meteorological data

The study area covers the province of Ghardaa (32.2° - 32.82° N, and 3.7° and 4.5° E), which is located in the desert region of Algeria (Figure 1), and at an altitude of 450m. Rainfall during the year is low in Ghardaia, it is classified as BWh climat [15]. In Ghardaia. The average precipitation is 68 mm, the average annual temperature is 21.0℃. The city has great solar potential throughout the year due to its location (the average daily solar radiation received is approximately 6000 Wh.m²) on a horizontal surface. The data sets used for this study include the total daily solar radiation on a horizontal surface measured and recorded at the Applied Research Unit for Renewable Energies URAER for the four-year period from January 1 (2013) to December 31 (2016). The first three years set of 2013-2015 was used as the training data set while the last year (2016) is used to test the different models.

2.2 Refinement of data

The accuracy of the models is greatly affected by the quality of the data used. It is preferable to perform the data cleaning procedure for improving the quality of the data by filtering it from any error or doubt.

The dataset of the daily solar radiation used includes unreliable values [16]. Therefore, we carried out a procedure in this work to filter the raw data before the design of the daily GSR.

1. For the knowledge of the inaccurate daily SR values, the daily clearness index K is calculated, the values which are outside the range 0.015 < K < 1 have been deleted [17].

2. A month is deleted from the dataset, if the incorrect values are greater than five days in this month; if the number is less than five, the values are replaced by correct values based on the interpolation [18]. Due to certain atmospheric phenomena such as cloud extinction and aerosol extinction that occur when solar radiation travels through the atmosphere, all values of H must be less than H0 in the data available, which means K < 1.

In Table 1, the variables used as inputs are very appropriate, where the close agreement is between Tmin, T_max, T_mean, H₀ and solar radiation. The lowest value of solar radiation is recorded during the month of December and the highest value is recorded in the month of July. In terms of relative humidity, this is the opposite of solar radiation, where the highest value is recorded during the month of July and the lowest value in the month of December.

The minimum temperature data and the minimum, medium and maximum humidity data are positively skewed with the average skewness factors of 0.02, 0.37, 0.24, 0.43 respectively, while the maximum and medium temperature data and extraterrestrial solar radiation, and daily incident solar radiation are negatively skewed with the average skewness factors of -0.04, -0.02, -0.24, -0.15, as expected.

1.png

Figure 1. The study area site

Table 1. The climatological cycle of daily global solar radiation for the study period

Statistically, two numerical measures of shape (skewness and excess kurtosis) can be used to test for normality. For kurtosis, the general guideline is that if the number is greater than +1, the distribution is peaked. If skewness is not close to zero, then your data set is not Gaussian distributed [19] as expected, all the data can be considered as Gaussian in their distributional behaviors. All the data can be considered as Gaussian in their distributional behaviors (Table 1).

2.3 Gaussian process regression (GPR)

GPR is a non-parametric model based on the Gaussian probability distribution [20]; it can be defined as a collection of random variables, of which any finite number GP has a joint Gaussian distribution [21]. Thus, a GP is completely specified by its 2^nd order statistics,

$f(x) \sim G P\left(m(x), k\left(x, x^{\prime}\right)\right)$    (1)

where, m(x) and $k\left(x, x^{\prime}\right)$ are the mean and covariance function of a real process f(x) respectively.

Suppose that a training set $\left\{\left(x_{i}, y_{i}\right), i=1 \ldots \ldots n\right\}$ The relationship between the $p-$ dimensional predictor $x \in \mathbb{R}^{p}$  and the target variable y is expressed as:

$y=f(x)+\varepsilon$   (2)

where, $\varepsilon$  is assumed to be an additive idd Gaussian noise, $\varepsilon \sim \mathrm{N}\left(0, \sigma_{n}^{2}\right)$.

The prior on the noisy observation becomes:

$\operatorname{cov}(y)=K(X, X)+\sigma_{n}^{2} I$    (3)

where, I denote the identy matrix of size n.

The joint distribution of the observed target values and the function values at the test locations prior is given by:

$\left|\begin{array}{l}

y \\

f_{*}

\end{array}\right| \sim \mathrm{N}\left(0 .\left|\begin{array}{cc}

K(X, X)+\sigma_{n}^{2} & K\left(X_{*}, X\right) \\

K\left(X_{*}, X\right) & K\left(X_{*}, X_{*}\right)

\end{array}\right|\right)$     (4)

where, $K\left(X, X_{*}\right)_{n \times n_{*}}$ denotes the covariance (or Gram) matrix between training test and also for different matrix $K\left(X_{*}, X\right), K\left(X_{*}, X_{*}\right) \operatorname{and} K(X, X) .$

The predictive equations for GPR becomes [22].

$f_{*} \mid X, y, X_{*} \sim \mathrm{N}\left(\bar{f}_{*}, \operatorname{cov}\left(f_{*}\right)\right)$    (5)

where,

$\operatorname{cov}\left(f_{*}=K\left(X_{*}, X_{*}\right)-\bar{f}_{*} K\left(X, X_{*}\right)\right.$    (6)

for a single test point $X_{*}$, the predictive distribution is a Gaussian distribution with mean and covariance given by:

$f_{*}=k_{*}^{T}\left(K+\sigma_{n}^{2} I\right)^{-1} y$     (7)

$\mathbb{V}\left[f_{*}\right]=k\left(x_{*}, x_{*}\right)-\bar{f}_{*} k_{*}$    (8)

where: $K=K(X, X), K_{*}=K\left(X, X_{*}\right) \text { and } k\left(x_{*}\right)=k_{*}$ denote the vector of covariances between the test point and the n training points. In the Eqns. (8) and (9), $\left(K+\sigma_{n}^{2} I\right)^{-1}$  can be calculated using Cholesky factorization [23].

2.4 Wavelet decomposition

The main motivation for using wavelet decomposition (WD) is the simple analysis of the series obtained. For many years, WD (or Wavelet transform) has been mixed with time series models as a preprocessing technique. WD uses a set of filters to decompose the original time series iteratively, so that separate forecasting models can be applied to each component. 

The continuous wavelet transform (CWT) of a function f(t), compared to the mother wavelet $\psi(t)$  can be written by the following integral [24]:

$F_{w}(a, \tau)=|a|^{-\frac{1}{2}} \int_{-\infty}^{+\infty} f(t) \psi^{*}\left(\frac{1-\tau}{a}\right) d t$    (9)

where, (*) represents the operation of the complex conjugation, $\tau \in \mathbb{R}$  is the translational value and $a \in \mathbb{R}^{+*}$  is the scaling coefficient. Unlike the Fourier transformation, the CWT has been discretized and is known as the discrete wavelet transform (DWT).

The approach is an implementation of the wavelet transform by scaling and translation of the wavelets in discrete time. In this case, the wavelets are given by:

$\psi_{n, k}(t)=\left|a_{0}^{n}\right|^{-\frac{1}{2}} \psi\left(\frac{1-k \tau_{0} a_{0}^{n}}{a_{0}^{n}}\right)$    (10)

where, n and k are integers and $a=a_{0}^{n}, \tau=k \tau_{0} a_{0}^{n}$.

More details on Wavelet transform can be found in the literature [24] and [25].

2.5 Structure of the hybrid model

We use W-GPR to predict daily solar radiation in the desert region. Through this study, we used wavelet analysis to decompose the time series of meteorological data into different components. The optimal GPR parameters are represented in the flow chart of Figure 2 based on the wavelet transform algorithm.

2.png

Figure 2. The flow chart of the proposed model of the wavelet-Gaussian process regression W-GPR

2.6 Model input data

Figure 3 shows the DWCs of the seven input variables using coiflet type wavelets (Table 2) with three levels of detailed decompositions and one level of approximation, where the approximate level has the lowest frequency. Many works use the entire wavelet sub-series [26], while other works delete the detail component and keep the remaining sub-series as noise based on the correlation coefficient [27].

In the proposed approach, we consider each wavelet-decomposed signal in its original form to capture their random attributes and their physical structure; on this basis, we insert the entire substring into the W-GPR model.

3.png

Figure 3. DWC of the inputs of the W-GPR model from 01-Jan-2013 to 31-December-2015 for Ghardaia Aero

Table 2. Effect of the wavelet type on model accuracy

2.7 Performance evaluation

The performance of the proposed models (W-GPR) are tested based on the following statistical measures:

$r^{2}=\frac{\left(\sum_{n=1}^{n}\left(H_{n, O b s}-\bar{H}_{n, \text { obs }}\right)\left(H_{n, \text { Pred }}-\bar{H}_{n, \text { Pred }}\right)\right)^{2}}{\sum_{i=1}^{N}\left(H_{n, O b s}-\bar{H}_{n, o b s}\right)^{2} \sum_{n=1}^{n}\left(H_{n, \text { Pred }}-\bar{H}_{n, \text { Pred }}\right)^{2}}$    (11)

$\mathrm{RMSE}=\sqrt{\frac{\sum_{\mathrm{n}=1}^{\mathrm{n}}\left(\mathrm{H}_{\mathrm{n}, \mathrm{Obs}}-\mathrm{H}_{\mathrm{n}, \mathrm{Pred})^{2}}\right.}{\mathrm{N}}}$    (12)

$\mathrm{MAE}=\frac{1}{\mathrm{~N}} \sum_{\mathrm{i}=1}^{\mathrm{N}}\left|\left(\mathrm{H}_{\mathrm{n}, \text { Pred }}-\mathrm{H}_{\mathrm{n}, \mathrm{Obs}}\right)\right|$    (13)

where,

H_i,Obs :observed values

H_n,Pred: predicted values

$\bar{H}$_n,obs: mean value of observations.

$\bar{H}$_n,pred: mean value of predictions.

N: total number of data.

r²: Coefficient of determination.

RMSE: Root Mean Square Error.

MAE: Mean Absolute Error.

rRMSE: relative Root Mean Square Error.

According to Paul et al. [28] the performance of the model by considering the rRMSE is defined as:

             rRMS E < 10 % the performance is Excellent.

10 % < rRMS E < 20 % the performance is Good.

20 % < rRMS E < 30 % the performance is Fair.

             rRMS E > 30 % the performance is Poor.

Since the prediction of solar radiation is very important in the management of solar systems. This paper investigated the possibility of prediction of W-GPR model, for daily solar radiation with high precision. Ten W-GPR models were developed using different combinations of inputs: Tmin, Tmax, RHmin, RHmax, RHmean and H₀. In order to evaluate the models and test their accuracy on the prediction, we used five statistical indicators. The results showed the significant effect of the wavelet type on the precision of the W-GPR models, where the wavelet type coif1 (Coiflet) has the highest precision, and the combination H₀, T_min, T_max, RH_min, RH_max, RH_mean offers great precision compared to the other proposed W-GPR models. To demonstrate the accuracy of the W-GPR model, its predictions are compared to the classical model (GPR). The results showed a significant improvement in the performance of the W-GPR model appearing in the statistical indices R² = 0.960, MAE=1.02 MJ/m²day, MSE=5.25 MJ/m²day, RMSE= 1.81 MJ/m²day, rRMSE= 11.21%. Finally, this model can be used to predict daily solar radiation in areas with a similar climate and can be further improved by introducing other variables, which should be the focus of our future work.