A Partially Amended Hybrid Bi-GRU—ARIMA Model (PAHM) for Predicting Solar Irradiance in Short and Very-Short Terms

Solar renewable energy (SRE) applications are substantial in eradicating the rising global energy shortages and reversing the approaching environmental apocalypse. Hence, effective solar irradiance forecasting models are crucial in utilizing SRE efficiently. This paper introduces a partially amended hybrid model (PAHM) by the implementation of a new algorithm. The algorithm innovatively utilizes bi-directional gated unit (Bi-GRU), autoregressive integrated moving average (ARIMA) and naive decomposition models to predict solar irradiance in 5-min and 60-min intervals. Meanwhile, the models’ generalizability strengths would be tested under an 11-fold cross-validation and are further classified according to their computational costs. The dataset consists of 32 months’ solar irradiance and weather conditions records. A fundamental result of this study was that the single models (Bi-GRU and ARIMA) outperformed the hybrid models (PAHM, classical hybrid model) in the 5-min predictions, negating the assumptions that hybrid models oust single models in every time interval. PAHM provided the highest accuracy level in the 60-min predictions and improved the accuracy levels of the classical hybrid model by 5%, on average. The single models were rigorous under the 11-fold cross-validation, performing well with different datasets; although the computational efficiency of the Bi-GRU model was, by far, the best among the models.


Significance and Problem Statement
Solar renewable energy (SRE) is the most abundant natural phenomenon occurring on Earth [1]. It has been deemed as the ultimate remedy for large scale environmental crises such as global warming, an effective alternative energy source for the depleting fossil fuels, and a viable business case for retaining profits. Therefore, SRE applications are common in residential, commercial, and industrial levels [2,3]. The international consensus is on devising flexible governmental policies for industries to invest in SRE and encouraging large scale deployments of solar complexes all over the world. As of 2019, more than 500 gigawatts of solar energy is produced worldwide; the investments are increasing almost by 100% yearly [4]. The Korean peninsula has a substantial photovoltaic power potential, resulting from the ample amount of incident global solar radiation, which is as high as 4.2 kWhm −2 [5]. Hence, there is significant economic feasibility and technical viability of SRE production in the region of interest, i.e., Jinju city, South Korea. On the other hand, the intermittent, inconsistent, and periodic nature of solar energy production counts as an inhibiting factor ahead. Seasonal changes of solar irradiance, as well as daily fluctuations in weather factors like the wind, temperature, and rain create highly unpredictable situations. This in turn complicates the integration of SRE into grids, hindering effective energy planning and scheduling [6].
Tremendous work has been carried out to tackle the uncertain characteristics of solar irradiance [7]. These studies have been applied to solar irradiance modeling, optimization, and short-or long-term predictions in order to integrate renewable energy smoothly into the grid. The developed models can be classified as probabilistic and statistical models, artificial intelligence (AI) applications, physical models, and hybrid models [8]. The hybridization approach combines two or three single methods such as AI and statistical hybrid models; benefiting from the strength of each model in building a final stronger model. In other words, time series datasets have inherent linear and non-linear characteristics. A decomposition method is used to decompose the time-series data to linear and non-linear segments [7]; consequently, single sub-models appropriate for linear or non-linear data are applied accordingly and individual results are combined to develop a final hybrid result. Thus, the merits of each model provide benefits in providing a final stronger result. So far, the hybrid models have outperformed the singles models in forecasting and modeling the solar irradiance [9]. Thus, in the current work, specifically the hybrid models are under focus. Additionally, a new approach of hybridization has been laid out in order to enhance accuracy levels in solar irradiance forecasts.

Related Literature
The authors in [10] have developed an effective hybrid model, consisting of wavelet decomposition (WD) and artificial neural network (ANN), for forecasting photovoltaic power. The high frequency of fluctuations in the time-series data were smoothed out by WD and then fed into the ANN model for final forecasts. This new hybrid method along with a clear classification of the sky clearness, has provided better RMSE values, between 0.093-0.229 (W), as compared to the single model. The ANN model with only three layers was used, which was shallow in structure, with only 5 days of input data under focus, was not large enough for training a machine learning algorithm thoroughly. Additionally, the developed algorithm could not predict fluctuations well as the inherent high frequency information in the time-series dataset was filtered out. In [11], an extensive study is delivered in utilizing seven machine learning (ML) models and developing three ensemble models for probabilistic forecasting of solar power generation. The dataset was grouped in hours and seven models were built based on each hour's data, which in turn, the corresponding forecasts were recalculated by the ensemble models; the latter models provided the best forecasts. Although the proposed method was very innovative but the broad amount of sheer calculations required, may incur unnecessary time losses, as similar predictions can be obtained with lesser calculations. The study in [8] provides a strong hybrid model composed of WD, convolutional neural network (CNN), and long short term memory model (LSTM). Firstly, a massive dataset was decomposed by WD and then the resulted smoothed and fluctuant parts of the data were used by to the CNN and LSTM models for further abstract feature extraction and temporal modeling, respectively. In general, the results of the hybrid model surpassed those of single models, i.e., the lowest average RMSE obtained was 100 (Wm −2 ). However, the errors in the forecasts were still high as compared to previous mentioned hybrid models; it also avoided in providing a generalized model for all of the seasons as a whole, as with a WD it was straightforward to build such a model. The authors of [12] presented a hybrid wind speed forecasting model comprised of WD, recurrent neural network (RNN) and adaptive neuro-fuzzy system (ANFIS). The model was new in applying ANFIS as the final step for providing step-wise wind predictions, instead of only averaging the final forecasts of the sub-models, as well as they had calculated probabilistic variances for determining accuracy intervals of the forecasts. The RMSE values obtained were as low as 0.968 (ms −1 ) for 15-min wind speed predictions. The model utilized RNN models only, which were fed with de-noised data; hence, some information were lost in this filtering process. This resulted in limiting the model in learning deep abstract features of the real time series dataset. In [13], a similar decomposition method as in [12] was applied along with neural networks and autoregressive models. Although this study was small, it reiterated the affectivity of the ensemble models. Overall, the smoothed and fluctuant components were better predicted by the combined models; the corresponding RMSE values were about 10.

Objectives
This paper proposes a partially amended hybrid model (PAHM) with the application of a new algorithm. The new algorithm combines the decomposition process, the deep RNN and ARIMA models, in an innovative way. Hence, the PAHM would be fundamentally different from the classical hybrid models. The deep RNN models which were developed in the previous work of this paper's author [14] would be utilized. In the mentioned work, a bi-directional gated recurrent unit (Bi-GRU) gave the best solar irradiance predictions; hence, for the first time the Bi-GRU model is proposed along with the ARIMA model, here. The naive decomposition (ND) will be applied to obtain fluctuant and smoothed time-series sub-datasets [7]. The Bi-GRU model would use the decomposed non-linear time-series dataset in addition to the actual data; while, the ARIMA model would use the linear decomposed time-series data. The ARIMA model will be developed by a trial and error process to select specific parameters for the single and hybrid cases, individually. The time series datasets under focus consist of weather conditions and solar irradiance. Additionally, the models' generalizabilities and strengths would be tested under a k-fold cross validation. The models would be classified further by comparing their computational efficiencies as well. In summary, the main objectives of this paper are: 1.
To develop a partially amended hybrid model (PAHM) consisting of deep RNN (Bi-GRU) with statistical ARIMA model.

2.
To develop an algorithm combining the decomposition process, deep RNN, and ARIMA models in a novel way.

3.
To improve the forecasting accuracies of solar irradiance in short terms as compared to the classical hybrid models.

4.
To test models with a k-fold cross validation and classify the models according to computational efficiencies.

Contents of the Paper
The rest of the paper is arranged as following: 2. Materials and Methods; 3. Results and Discussion; 4. Conclusions. The last parts are comprised of Acknowledgments and References.

Deep RNN Architecture Development
Recurrent neural networks (RNN) are widely applied in non-linear time series problems, speech recognition, sentiment classification, modeling, or predictions [15]. They have the ability to memorize data, therefore, are suitable in learning the past characteristics of a time series dataset and utilize it in future sequence analysis [16]. A simple RNN architecture is shown in Figure 1a,b; basically, it is a recursive model which predicts h t , state at time t, using information from past states h t−1 with a differentiable function f. The corresponding weights, U, V, and W, are shared across the model and are governed by Equations (1)-(3) [17].
On the other hand, the gated recurrent unit (GRU) is a developed variant of RNN models. It has a complex structure as compared to its predecessors and copes better against issues persistent in classical Energies 2020, 13, 435 4 of 20 RNN models such as vanishing gradients [18]. A GRU cell has internal gates with write, read, and reset capabilities, i.e., it can store and update information by manipulating these gates. A GRU cell is depicted in Figure 1d and the corresponding equations are as following [19]. In the equations, b, g, R, and σ denotes the bias vectors, non-linear activation functions, weight parameters and sigmoid functions, respectively.
reset gate : candidate state : new state : As explained thoroughly in the previous work of this paper's author [14], bidirectional RNN (Bi-RNN) models are very effective in predicting short term solar irradiance. These models can carry knowledge from the past sequences as well as future sequences of a time series, hence they are called bidirectional. The Bi-RNN models have forward and backward recurrent hidden layers with same inputs and outputs. In a time series of sequence T, the forward recurrent layer calculates hidden states from t = 1 to t = T; while the backward recurrent layer calculates hidden states from t = T to t = 1. Hence, in Bi-RNN models the entire fundamental RNN processes are performed twice, which consequently increases the computational costs; but the affectivity of Bi-RNN models deem them more preferable than classical unidirectional ones [19]. Lastly, the results of the forward and backward predictions are merged to get a final prediction. A simple architecture of a Bi-RNN model is depicted in Figure 1d, and the underlying relationships are demonstrated in Equations (9) and (10) [19].
Energies 2020, 13, x FOR PEER REVIEW 4 of 21 classical RNN models such as vanishing gradients [18]. A GRU cell has internal gates with write, read, and reset capabilities, i.e., it can store and update information by manipulating these gates. A GRU cell is depicted in Figure 1d and the corresponding equations are as following [19]. In the equations, b, g, R, and denotes the bias vectors, non-linear activation functions, weight parameters and sigmoid functions, respectively.
: ℎ = ℎ − 1 * : : ℎ = (1 − ) * ℎ − 1 + * As explained thoroughly in the previous work of this paper's author [14], bidirectional RNN (Bi-RNN) models are very effective in predicting short term solar irradiance. These models can carry knowledge from the past sequences as well as future sequences of a time series, hence they are called bidirectional. The Bi-RNN models have forward and backward recurrent hidden layers with same inputs and outputs. In a time series of sequence T, the forward recurrent layer calculates hidden states from = 1 to = ; while the backward recurrent layer calculates hidden states from = to = 1. Hence, in Bi-RNN models the entire fundamental RNN processes are performed twice, which consequently increases the computational costs; but the affectivity of Bi-RNN models deem them more preferable than classical unidirectional ones [19]. Lastly, the results of the forward and backward predictions are merged to get a final prediction. A simple architecture of a Bi-RNN model is depicted in Figure 1d, and the underlying relationships are demonstrated in Equations (9) and (10) [19]. As mentioned earlier, the Bi-RNN models exploit the availability of training data sequences as a whole; while learning about a sequence at time t, it considers the hidden states prior and post that sequence. Since weather time series sequences are highly correlated to the immediate past and future sequences, the model utilizes this information in learning about the abstract features of the dataset [20]. It is also pointed out in the literature that datasets with high fluctuations and less causal relationship are harder for unidirectional RNN models to predict [21]. There are two classes of Bi-RNN models which are bidirectional gated recurrent networks (Bi-GRU) and bidirectional long-short term memory (Bi-LSTM). In this study, the former type is selected, as its superior performance, in predicting solar irradiance, over the latter type has been shown in [14]. Table 1 shows the respective hyper-parameters like transfer or activation function, epoch, optimizer, number of neurons, and hidden layers were selected after trial and error processes in [14]. However, in the current work, the architectural parameters of the Bi-GRU model would be changed as per the new datasets in order to provide reasonable accuracy levels. As mentioned earlier, the Bi-RNN models exploit the availability of training data sequences as a whole; while learning about a sequence at time t, it considers the hidden states prior and post that sequence. Since weather time series sequences are highly correlated to the immediate past and future sequences, the model utilizes this information in learning about the abstract features of the dataset [20]. It is also pointed out in the literature that datasets with high fluctuations and less causal relationship are harder for unidirectional RNN models to predict [21]. There are two classes of Bi-RNN models which are bidirectional gated recurrent networks (Bi-GRU) and bidirectional long-short term memory (Bi-LSTM). In this study, the former type is selected, as its superior performance, in predicting solar irradiance, over the latter type has been shown in [14]. Table 1 shows the respective hyper-parameters like transfer or activation function, epoch, optimizer, number of neurons, and hidden layers were selected after trial and error processes in [14]. However, in the current work, the architectural parameters of the Bi-GRU model would be changed as per the new datasets in order to provide reasonable accuracy levels.

ARIMA Model Description
Autoregressive integrated moving average (ARIMA) is an effective model in non-stationary time series forecasting [22]. It is influential in building on the past characteristics of a time series and predict its future direction [23]. In an ARIMA model the predictions of a time series, y t , are assumed to be a linear function of the past observations and random errors, ε t . Generally, an ARIMA model can be expressed as [19]  Here p, d, and q are integers indicating the lag order of autoregressive (AR), degree of differencing (I) and size of the moving average (MA) parts of the ARIMA model, respectively, and it is denoted as ARIMA (p,d,q). On the other hand, L shows the lag or backshift operator and ϕ(L) is the lag polynomial of the model. In order to apply the ARIMA model, first the non-stationary time series is changed to stationary time series by finite differencing of the observations. The differencing operation would remove instability in the time series and result in constant statistical mean and autocorrelations [24,25].
In the literature [26], the prediction process of the ARIMA model is divided into three steps: identification, estimation, and evaluation; where d is selected in the first step and the values of p and q are estimated in the second step. The value of d is selected 1 in most cases [27] and in this study, the mentioned value would be applied. In the estimation step, the autocorrelation function (ACF) and partial autocorrelation function (PACF) are proposed to determine p and q orders of the ARIMA model in order to create stationarity. By applying different orders of p and q, the ACF and PACF are plotted against each other in order to check the level of irregularity-i.e., existing trend-in the time series. Consequently, the orders of p and q are estimated from the plots having no or ignorable level of trend [23]. Three performance criteria would be utilized in order to evaluate the performance of the ARIMA model which is discussed later in this section. The development of the ARIMA model is usually followed by trial and error of the above steps until a plausible model is estimated.

Naive Decomposition Method
Generally, non-stationary time series demonstrates certain characteristics of trend, cycle, seasonality, and irregularity. Hence, decomposition methods are widely used for splitting a non-stationary time series dataset into constituent parts, i.e., smoothed and fluctuant parts [28]. Consequently, each segment is analyzed by a suitable well-established method. In this study, the naive decomposition method is used to obtain stationary and non-stationary datasets from the solar irradiance and weather time series. There are other decomposition methods like wavelet decomposition [29] or empirical mode decompositions [30], but these models are basically derived from the naive classical decomposition. Therefore, in this study, the naive method is sufficient.
Naive decomposition is also called 'classical decomposition' due to its simplicity [31]. This method splits a time series into trend (smoothed) and residual (fluctuant or noise) components [32]. There are additive and multiplicative types of naive decomposition; the former type is applicable on time series where seasonal variations over the years are almost constant [33]. In this study, the additive model will be applied as a 3-year dataset of solar irradiance and demonstrates negligible year by year variations, i.e., the seasonal changes between same seasons of different years will be considered zero. The classical decomposition method splits a time series, y t , into trend, L t , (linear or smoothed) and residual, R t , (fluctuant or noise) components.

Data Preparation
The time series dataset would consist records of solar irradiance (Wm −2 ), sun hour, temperature ( • C), relative humidity (%), wind speed (ms −1 ), and direction ( • ). The sequences will be recorded and averaged over 5 min and 1 h periods and the data between 05:00-19:00 will be included in the analysis. The corresponding weather station, shown in Figure 2, is located in the Gyeongsang National University, JinjuCity, Republic of South Korea. It is equipped with sunshine sensor, pyrometer collecting solar radiation, weather sensors, and solar panel for obtainable solar energy measurements.

Data Preparation
The time series dataset would consist records of solar irradiance (Wm −2 ), sun hour, temperature (°C), relative humidity (%), wind speed (ms −1 ), and direction (°). The sequences will be recorded and averaged over 5 min and 1 h periods and the data between 05:00-19:00 will be included in the analysis. The corresponding weather station, shown in Figure 2, is located in the Gyeongsang National University, JinjuCity, Republic of South Korea. It is equipped with sunshine sensor, pyrometer collecting solar radiation, weather sensors, and solar panel for obtainable solar energy measurements. The time series data will be preprocessed before applying to the single and hybrid models. Raw data may contain outliers, missing values, or inconsistencies. The outliers and missing values would be changed with an average value of the immediate post and prior sequences or would be deleted completely. Some inconsistencies may stem from sensor abnormalities. Therefore, percentile distributions, standard deviations, and means of each components in the data would be examined to verify the dataset's logicality and sensibility.

Classical Hybrid Model
In the classical hybrid models, the time series dataset would be decomposed to linear, , and non-linear, , parts. Accordingly, each component would be used by one sub-model to provide component-wise predictions [28]. The final predictions of each sub-model would be averaged to get the final solar irradiance predictions and consequently the losses , as shown in Equations (14) The time series data will be preprocessed before applying to the single and hybrid models. Raw data may contain outliers, missing values, or inconsistencies. The outliers and missing values would be changed with an average value of the immediate post and prior sequences or would be deleted completely. Some inconsistencies may stem from sensor abnormalities. Therefore, percentile distributions, standard deviations, and means of each components in the data would be examined to verify the dataset's logicality and sensibility.

Classical Hybrid Model
In the classical hybrid models, the time series dataset would be decomposed to linear, L t , and non-linear, R t , parts. Accordingly, each component would be used by one sub-model to provide component-wise predictions [28]. The final predictions of each sub-model would be averaged to get the final solar irradiance predictionsŷ t and consequently the losses e t , as shown in Equations (14) and (15). In this study, the former component would be used by the ARIMA model and the latter component would be fed into the Bi-GRU model for solar irradiance predictions. This is the classical method followed in the literature [22,25,26]. This method would provide a basis for comparison with the new hybrid model proposed in the current work.

Partially Amended Hybrid Model (PAHM)
Here, a partially amended hybrid model (PAHM) would be developed by using the following new algorithm. This algorithm combines the hybrid model in a novel way. While in classical hybrid models decomposition results are directly utilized by sub-models [28], in this approach the actual data, y t , and the components of the decomposed data will be used in a different way by the sub-models. The decomposition method provides smoothed, L t , and fluctuant, R t , new time series datasets from the actual series. As mentioned before, the weather conditions in the Jinju city has high irregularities. This makes modeling and forecasting hard for RNN models. Hence, R t would provide an effective information of the weather uncertainties to the Bi-GRU sub-model. While the ARIMA sub-model would use the smoothed dataset, the deep RNN model would be fed with the actual time series as well as the fluctuant residual data set from the decomposition. The final prediction,ŷ t , would be obtained from by averaging the outputs of the two sub-models. In Figure 3, the new proposed algorithm flowchart is depicted. and (15). In this study, the former component would be used by the ARIMA model and the latter component would be fed into the Bi-GRU model for solar irradiance predictions. This is the classical method followed in the literature [22,25,26]. This method would provide a basis for comparison with the new hybrid model proposed in the current work.

Partially Amended Hybrid Model (PAHM)
Here, a partially amended hybrid model (PAHM) would be developed by using the following new algorithm. This algorithm combines the hybrid model in a novel way. While in classical hybrid models decomposition results are directly utilized by sub-models [28], in this approach the actual data, yt, and the components of the decomposed data will be used in a different way by the submodels. The decomposition method provides smoothed, Lt, and fluctuant, Rt, new time series datasets from the actual series. As mentioned before, the weather conditions in the Jinju city has high irregularities. This makes modeling and forecasting hard for RNN models. Hence, Rt would provide an effective information of the weather uncertainties to the Bi-GRU sub-model. While the ARIMA sub-model would use the smoothed dataset, the deep RNN model would be fed with the actual time series as well as the fluctuant residual data set from the decomposition. The final prediction, , would be obtained from by averaging the outputs of the two sub-models. In Figure 3  The root mean square error (RMSE), coefficient of determination (R 2 ), and mean absolute error (MAE) metrics will be used to demonstrate the accuracy level of the models. These metrics are important in the development of the new hybrid model as the time series dataset under consideration possesses highly volatile data. The RMSE performance metric is sensitive to bigger errors in the predictions and the R 2 criterion depicts how fit a model is to the actual non-linear time series; meanwhile, the MAE level would shed light on the average distribution of errors in the model predictions as a whole. The first metric is a good indicator of the strength of the model against high perturbations, while the second one demonstrates the level of bias in the predictions and the last one shows how much the forecasted values are spread over the entire model [34]. The formulas of these metrics are as below. The inherent weather conditions demonstrate complex characteristics and high fluctuations from season to season [35]. For example, temperature changes between −13 • C to 38 • C and a standard deviation of 10 • C, during a year, as shown in Table 2. Hence, Bi-GRU and ARIMA models' generalizability needs to be scrutinized; an 11-fold cross validation would be utilized to check each single model's performance in each season. Additionally, this method would solve the under-fitting and over-fitting issues of the mentioned models [36]. Computational Costs of the Models There will be different time requirements while training and fitting the Bi-GRU and ARIMA models. In case of deep learning models, as the number of neurons and hidden layers increase, so does the required time for training these models [24]. On the other hand, the statistical based ARIMA model would consume more time as the time series dataset gets bigger. Consequently, the time efficiencies of each single model would be studied. It would be critical to go beyond just comparing the accuracy levels of the models as computational efficiencies provides an effective criterion while selecting the best model for predicting the solar irradiance. The computational efficiency would be inferred from the amount of time consumed for training and testing a model. In this work, all the models would be built in Python and run on a computer with the following specifications.

Time Series Dataset Preprocessing
The solar irradiance and weather factors were collected over a period of 32 months, starting from March 2017. The sequences were reduced to the daily sun hour intervals i.e., from 05:00 to 19:00. In the resulted time series dataset, there were 174,914 and 14,756 records for the 5-min and 1 h periods, respectively.

Data Checkup
In the reshaped time series datasets, there were six outlier values which were seemingly unfitting the corresponding trends of those particular instances; while, two missing values in the sequences were identified. The unfitting and missing instances were changed to the average of the values in the immediate vicinity of those particular records. The corresponding amendment was minimal and would not affect the final calculations. Additionally, the general logicality of the datasets were depicted by the statistical measures shown in Table 2a,b; the percentile distributions are evenly increasing from minimum to maximum. Meanwhile, the standard deviations and the means for all variables demonstrated higher accumulation of data points around the mean [37]. After finishing the preprocessing, the time series datasets for the two periods were ready for use in the models.

Variable Characteristics
The correlation analysis of the variables in the time series datasets exhibited the levels of sensitivity between solar irradiance and other weather factors. The air temperature, wind speed, relative humidity, and sun hour variables had higher correlation values. In the 5-min time series dataset, the temperature correlation value is below 0.5. Nevertheless, all of these variables were used, as it is deemed reasonable to feed the neural networks with as much higher amounts of data as possible [38]. The corresponding correlation levels in the 5-min and 60-min periods, are depicted in Table 3.

Single Model Results
It was important to select proper architectures for the ARIMA and Bi-GRU models in order to obtain reliable results. After an extensive trial and error processes, specific model parameters were selected. Consequently, the two models were applied on the processed time series datasets for very short term (5-min) and short term (60-min) solar irradiance predictions.
The Bi-GRU model specifications given in Table 1, were amended to fit the new 32 months 5-min dataset; the number of hidden layers and total neurons were increased, as it was necessary for the deep RNN model to learn proportionally new abstract features of the bigger datasets [14]. As for the 60 min solar irradiance predictions, the corresponding dataset was smaller; hence, a shallower RNN architecture was selected, i.e., fewer neurons and hidden layers. The respective configurations of the Bi-GRU model for 5-min and 60-min intervals were chosen as shown in Table 4. On the other hand, for the ARIMA model, the lag order, degree of difference and the size of the moving average were selected 3, 1, and 0, respectively. The value 1 was selected for the degree of differencing, as it is appropriate for most time series datasets. The corresponding ACF and PACF plots, in Figure 4c,d, showed sinusoidal characteristics and approached zero faster after lag 1. Therefore, the time series depicted AR properties with possible orders between 3 and 4 and an MA value of 0. After running the ARIMA model for both lag orders, 3 was selected; the mean Gaussian distribution of errors for both time intervals, Figure 4a,b, were nearest to zero with the mentioned p, d, q values. The distribution of errors concentrating around zero, demonstrated that the residuals after the differencing process did not possess any trend; hence, the transformed time series showed the highest degree of stationarity with the mentioned p, d, and q. This is a required condition in showing the stationarity of the time series and building any ARIMA model [21].

Single Model Results
It was important to select proper architectures for the ARIMA and Bi-GRU models in order to obtain reliable results. After an extensive trial and error processes, specific model parameters were selected. Consequently, the two models were applied on the processed time series datasets for very short term (5-min) and short term (60-min) solar irradiance predictions.
The Bi-GRU model specifications given in Table 1, were amended to fit the new 32 months 5min dataset; the number of hidden layers and total neurons were increased, as it was necessary for the deep RNN model to learn proportionally new abstract features of the bigger datasets [14]. As for the 60 min solar irradiance predictions, the corresponding dataset was smaller; hence, a shallower RNN architecture was selected, i.e., fewer neurons and hidden layers. The respective configurations of the Bi-GRU model for 5-min and 60-min intervals were chosen as shown in Table 4. On the other hand, for the ARIMA model, the lag order, degree of difference and the size of the moving average were selected 3, 1, and 0, respectively. The value 1 was selected for the degree of differencing, as it is appropriate for most time series datasets. The corresponding ACF and PACF plots, in Figure 4c,d, showed sinusoidal characteristics and approached zero faster after lag 1. Therefore, the time series depicted AR properties with possible orders between 3 and 4 and an MA value of 0. After running the ARIMA model for both lag orders, 3 was selected; the mean Gaussian distribution of errors for both time intervals, Figure 4a,b, were nearest to zero with the mentioned p, d, q values. The distribution of errors concentrating around zero, demonstrated that the residuals after the differencing process did not possess any trend; hence, the transformed time series showed the highest degree of stationarity with the mentioned p, d, and q. This is a required condition in showing the stationarity of the time series and building any ARIMA model [21]. The respective RMSE, R 2 , and MAE performance results of the models are presented in Table 5. The single models performed similar in the mentioned intervals, while they were better in the very short term solar irradiance predictions. This was partly due to the fact that the corresponding dataset for the very short term was very large, which resulted in effective training of the models. Additionally, the inherent characteristics of the datasets, depicted in Figure 5a,b, vary for the individual time intervals. The solar irradiance has high variability-i.e., volatility-between consecutive values in the 60-min dataset as compared to the 5-min ones. Hence, forecasting one step ahead is harder for the models in the longer term, i.e., 60 min.
These results demonstrate the limitations of the Bi-GRU and ARIMA models as they behave better only with huge datasets which possess a higher degree of stationarity [39]. However, in general, ARIMA is a proven method for forecasting non-stationary time series and solves many limitations of its predecessors [26]; while GRU is a state of the art neural network with internal The respective RMSE, R 2 , and MAE performance results of the models are presented in Table 5. The single models performed similar in the mentioned intervals, while they were better in the very short term solar irradiance predictions. This was partly due to the fact that the corresponding dataset for the very short term was very large, which resulted in effective training of the models. Additionally, the inherent characteristics of the datasets, depicted in Figure 5a,b, vary for the individual time intervals. The solar irradiance has high variability-i.e., volatility-between consecutive values in the 60-min dataset as compared to the 5-min ones. Hence, forecasting one step ahead is harder for the models in the longer term, i.e., 60 min.
These results demonstrate the limitations of the Bi-GRU and ARIMA models as they behave better only with huge datasets which possess a higher degree of stationarity [39]. However, in general, ARIMA is a proven method for forecasting non-stationary time series and solves many limitations of its predecessors [26]; while GRU is a state of the art neural network with internal memory cells which can learn and maintain very abstract features of any time series such solar irradiance temporal series [24]. As an example, the predictions of the ARIMA and Bi-GRU models were compared with the actual solar irradiance values and depicted in Figure 5a,b. As seen in the graphs, the performance levels of the models were similar in both of the intervals of the solar irradiance forecasts. During the selected days, the actual values of solar irradiance were fluctuating violently for both short and very short term intervals. Nonetheless, both of the models were rigorous in predicting solar irradiance levels of the highly turbulent days, as well as in the quite day. However, the RMSE and MAE values were still high due the fluctuating characteristics of the solar irradiance time series dataset for both short (60-min) and very short (5-min) terms.
Energies 2020, 13, x FOR PEER REVIEW 13 of 21 memory cells which can learn and maintain very abstract features of any time series such solar irradiance temporal series [24]. As an example, the predictions of the ARIMA and Bi-GRU models were compared with the actual solar irradiance values and depicted in Figure 5a,b. As seen in the graphs, the performance levels of the models were similar in both of the intervals of the solar irradiance forecasts. During the selected days, the actual values of solar irradiance were fluctuating violently for both short and very short term intervals. Nonetheless, both of the models were rigorous in predicting solar irradiance levels of the highly turbulent days, as well as in the quite day. However, the RMSE and MAE values were still high due the fluctuating characteristics of the solar irradiance time series dataset for both short (60-min) and very short (5-min) terms.

K-Fold Cross Validation Results of the Single Models
There were 11 seasons in the time series datasets, so an 11-fold cross validation was applied on the ARIMA and Bi-GRU models. The K-fold cross validation results of the models were similar; while both models performed better in very short term solar irradiance predictions. The average RMSE values of the Bi-GRU model over the 11 seasons were 61.64 and 100.14 W/m 2 for 5-min and 60-min solar irradiance predictions, respectively. While those of the ARIMA model were 59.43 and 97.52 W/m 2 for 5 min and 60 min solar irradiance predictions, respectively. The seasonal average RMSE, R 2 , and MAE values of each model are given in the scatter plots in Figure 6a-c.
Overall, the prediction accuracies of the models differed for each season; both models performed better in the winter while provided the worst accuracy levels in the summer seasons. This indicates the fluctuating nature of the weather factors as well as the solar irradiance levels throughout the seasons. The higher predictability of the winter seasons might stem from the sun hours difference. In the winter, the sun hours are approximately 4 h less than those of the summer seasons. Hence, it was harder for the models to learn and predict the solar irradiance in the summer season.

Computational Cost Efficiencies
The time consumed by each single model was collected to compare their time efficiencies. As seen in Table 6, the Bi-GRU model required a very low amount of time in the training and testing phases, as compared to the ARIMA model. Hence, the time efficiency of the former model has put it in an advantageous position against the counterpart model. On the other hand, the hybrid models would constitute both of these models. Therefore, their time efficiencies would be even higher than individual single models, as a whole.
Energies 2020, 13, x FOR PEER REVIEW 14 of 21 There were 11 seasons in the time series datasets, so an 11-fold cross validation was applied on the ARIMA and Bi-GRU models. The K-fold cross validation results of the models were similar; while both models performed better in very short term solar irradiance predictions. The average RMSE values of the Bi-GRU model over the 11 seasons were 61.64 and 100.14 W/m 2 for 5-min and 60-min solar irradiance predictions, respectively. While those of the ARIMA model were 59.43 and 97.52 W/m 2 for 5 min and 60 min solar irradiance predictions, respectively. The seasonal average RMSE, R 2 , and MAE values of each model are given in the scatter plots in Figure 6a-c.
Overall, the prediction accuracies of the models differed for each season; both models performed better in the winter while provided the worst accuracy levels in the summer seasons. This indicates the fluctuating nature of the weather factors as well as the solar irradiance levels throughout the seasons. The higher predictability of the winter seasons might stem from the sun hours difference. In the winter, the sun hours are approximately 4 h less than those of the summer seasons. Hence, it was harder for the models to learn and predict the solar irradiance in the summer season.

Computational Cost Efficiencies
The time consumed by each single model was collected to compare their time efficiencies. As seen in Table 6, the Bi-GRU model required a very low amount of time in the training and testing phases, as compared to the ARIMA model. Hence, the time efficiency of the former model has put it in an advantageous position against the counterpart model. On the other hand, the hybrid models would constitute both of these models. Therefore, their time efficiencies would be even higher than individual single models, as a whole.

Classical Hybrid Model Performance
Firstly, the hybrid model was applied in the classical way, i.e., without the proposed algorithm. The non-stationary time series dataset was decomposed to stationary and non-stationary components by a naive additional decomposition. As seen in the graphs in Figure 7a Table 7, the classical hybrid model performed better in the 60-min solar irradiance predictions while it provided higher errors in the 5-min solar predictions.

Classical Hybrid Model Performance
Firstly, the hybrid model was applied in the classical way, i.e., without the proposed algorithm. The non-stationary time series dataset was decomposed to stationary and non-stationary components by a naive additional decomposition. As seen in the graphs in Figure 7a Table 7, the classical hybrid model performed better in the 60-min solar irradiance predictions while it provided higher errors in the 5-min solar predictions.

PAHM Performances
In this part, the solar irradiance prediction results of the new hybrid model, amended with the new proposed algorithm, would be demonstrated. As stated in Section 2.5, the time series dataset was decomposed, as shown in Figure 7a Table 8.

PAHM Performances
In this part, the solar irradiance prediction results of the new hybrid model, amended with the new proposed algorithm, would be demonstrated. As stated in Section 2.5, the time series dataset was decomposed, as shown in Figure 7a

Final Comparisons of Different Models
So far, the single and hybrid models provided mixed results. Table 9, demonstrate the level of improvements by the hybrid models as compared to the single models. There were significant forecast accuracy improvements by the hybrid models in the 60-min solar irradiance predictions. While the hybrid models provided less satisfactory results in the 5-min time step, as compared to the single models.
These results might stem from the fact that the 5-min solar irradiance time series dataset has higher standard deviations in the weather factors, as depicted in Table 2. Meanwhile the 60-min dataset looks more stable in this respect. Additionally, the former dataset does not have the Sun hour data, while it is included in the latter dataset; hence, predicting solar irradiance would be easier for the hybrid model in this time step. It can be deduced that hybrid models would not be very effective in predicting solar irradiance in very short time steps, i.e., 1-to 5-min time steps.
On the other hand, the inclusion of the new algorithm to the hybrid model has improved the forecasting accuracies. As seen in Table 9, the new hybrid model, with the new algorithm, performed better than classical hybrid model in both time steps. It can be concluded that the addition of an actual time series dataset along with the decomposed non-linearity information enhances the prediction abilities of the Bi-GRU model, and hence that of the hybrid model.In Figure 8, the classical hybrid model and the PAHM results are compared in the two time intervals. Both hybrid models performed better in the 60-min short-term predictions, while the 5-min very short-term predictions are far from being satisfactory and lower than those of the single models.

Final Comparisons of Different Models
So far, the single and hybrid models provided mixed results. Table 9, demonstrate the level of improvements by the hybrid models as compared to the single models. There were significant forecast accuracy improvements by the hybrid models in the 60-min solar irradiance predictions. While the hybrid models provided less satisfactory results in the 5-min time step, as compared to the single models.
These results might stem from the fact that the 5-min solar irradiance time series dataset has higher standard deviations in the weather factors, as depicted in Table 2. Meanwhile the 60-min dataset looks more stable in this respect. Additionally, the former dataset does not have the Sun hour data, while it is included in the latter dataset; hence, predicting solar irradiance would be easier for the hybrid model in this time step. It can be deduced that hybrid models would not be very effective in predicting solar irradiance in very short time steps, i.e., 1-to 5-min time steps.
On the other hand, the inclusion of the new algorithm to the hybrid model has improved the forecasting accuracies. As seen in Table 9, the new hybrid model, with the new algorithm, performed better than classical hybrid model in both time steps. It can be concluded that the addition of an actual time series dataset along with the decomposed non-linearity information enhances the prediction abilities of the Bi-GRU model, and hence that of the hybrid model.In Figure 8, the classical hybrid model and the PAHM results are compared in the two time intervals. Both hybrid models performed better in the 60-min short-term predictions, while the 5-min very short-term predictions are far from being satisfactory and lower than those of the single models. It can be seen from these results that the single models are preferable in the very short term solar irradiance predictions. While the ARIMA and Bi-GRU models provided similar accuracy levels in the 5-min predictions, the Bi-GRU model can be deemed more suitable as it had higher computational efficiency too. On the other hand, the hybrid model with new algorithm provided higher accuracy levels in the longer time steps, i.e., 60-min. Table 9. Performance improvements by the hybrid models as compared to the single models. It can be seen from these results that the single models are preferable in the very short term solar irradiance predictions. While the ARIMA and Bi-GRU models provided similar accuracy levels in the 5-min predictions, the Bi-GRU model can be deemed more suitable as it had higher computational Energies 2020, 13, 435 18 of 20 efficiency too. On the other hand, the hybrid model with new algorithm provided higher accuracy levels in the longer time steps, i.e., 60-min.

Conclusions
Solar renewable energy (SRE) is considered as a solution in the face of the perceived environmental issues, finishing fossil fuels, and the ever-increasing energy demand of the industries. Hence, it is important to develop models that would significantly ease the integration of the SRE. This paper developed a partially amended hybrid model (PAHM) by the introduction of a new algorithm which combined the Bi-GRU, ARIMA, and the naive decomposition processes in an innovative way. The single models (Bi-GRU, ARIMA) and the hybrid models (PAHM, classical hybrid model) were applied over a 32-month time series dataset, to forecast solar irradiance in 5-min and 60-min intervals. These models were further classified based on an 11-fold cross validation and computational cost efficiency levels. The single models' results (with an average RMSE of 71.17 W/m 2 ) were better in the 5-min solar irradiance predictions while the PAHM outperformed other models in the 60-min time step (with an RMSE of 52.64 W/m 2 ). The former results were contrary to hypothesis of the hybrid models outperforming single models at any time interval. Furthermore, the proposed new algorithm was effective in enhancing the accuracy levels of the classical hybrid model, on average by 5%, in both time intervals. On the other hand, under the 11-fold cross validation, ARIMA and Bi-GRU models performed similarly, showing that both models were generalizable and robust with different datasets. However, there were huge differences between the Bi-GRU and other models in terms of computational time efficiencies. The Bi-GRU model consumed significantly less time in the training and testing phase on the same dataset, as compared to the ARIMA and hybrid models. This study can be established further by using weather data from different geographical locations. Moreover, the utilization of different decomposition methods can enhance the generalizability of the findings of PAHM. In the future, these points would be considered to promote the valuable findings of this work.