1. Introduction

Water resources’s overexploitation associated with global warming affects majorly the quality and amount of water used (^{Yurekli, 2021}). Agricultural water demand is supplied by planning water resources and their correct use. Therefore, evapotranspiration (ET) needs should be determined for sustainable water management. ET is the most important factor that represents the effect of climate parameters on the hydrologic cycle (^{Pereira et al., 2015}). Direct calculation methods of ET are not practical nor economical. Instead, reference evapotranspiration (ET₀) is successfully estimated using climate parameters with the development of empirical equations. The Penman-Monteith method is proposed by the Food & Agriculture Organization (FAO) as the best-fit method for calculating ET₀. The changes in hydrological and meteorological components affect ET₀ seriously. The most significant problem in calculating ET₀ is the missing data. ^{Allen et al. (1998)} reported that instead of calculating alternative ET₀ with limited data, it is more appropriate to calculate ET₀ by estimating the missing data.

Meteorology stations in many regions of the world do not have enough data to make ET₀ calculations with the Penman-Monteith method (^{Zanetti et al., 2019}). Solar radiation data representing energy from the sun on ET₀ is one of the most commonly lacking meteorological parameters. Also, solar radiation is a significant variable since it is used in various areas such as meteorology, solar energy systems (^{Yadav and Chandel, 2014}), and hydrological and agricultural studies (^{Yang et al., 2006}). The missing solar radiation data can be estimated from the sunshine duration using the ^{Angstrom equation (Angstrom, 1924)}, and for regions where the sunshine duration is not measured, from the maximum and minimum air temperature difference reported by ^{Hargreaves and Samani (1982)}. Artificial intelligence techniques (^{Jiang, 2009}; ^{Yacef et al., 2012}) and machine learning models (^{Lauret et al., 2015}) used in solar radiation estimation give more successful results than empirical equations. While some researchers used the empirical equations for solar radiation estimates (^{Besharat et al., 2013}; ^{Hassan et al., 2016}; ^{Yıldırım et al., 2018}, ^{Mohammadi and Moazenzadeh, 2021}), some researchers employed machine learning methods based on the relationship of inputs and outputs learned from previous data sets (^{Sözen et al., 2004}; ^{Amrouche and Le Pivert 2014}; ^{Shamshirband et al., 2015}; ^{Premalatha and Valan Arasu 2016}; ^{Guermoui and Rabehi 2018}).

^{Sözen et al. (2005)} reported higher success rates with the ANN method compared to the classic regression models for solar radiation estimation in Turkey. ^{Rahimikhoob (2010)}, in his study conducted in Iran, modeled solar radiation using the Hargreaves equation and the ANN method obtaining the best estimate with ANN. Compared to the empirical equations, the SVM model gave better results in places where air temperature data are available (^{Chen et al., 2011}). ^{Premalatha and Valan Arasu (2016)} stated that Levenberg-Marquardt was the ANN algorithm that gave the best result for solar radiation estimation in India. On the other hand, ^{Chandola et al. (2020)} compared LSTM and RNN methods for solar radiation estimation in the Thar desert of Pakistan and obtained the most successful results from the LSTM model.

The purpose of the present study was to estimate and analyze, using deep learning and machine learning models, the solar radiation values measured for Antalya, Mersin, Mugla stations, which receive the highest radiation levels in Turkey, and Rize, Trabzon, and Ordu stations, which receive the lowest radiation levels. In addition, the performances of empirical models such as Hargreaves-Samani and Angstrom equations were compared with LSTM, SVM, and ANN models. In the present study, LSTM, ANN, and SVM algorithms were developed using the inputs of two different empirical models. The ability to perform estimations with unlimited variables by learning historical data enables these models to be used frequently. The suitability for modeling dynamic environments such as climate parameters, minimizing human involvement, high learning speed, and accuracy made learning models preferable to empirical models.

2. Material and Method

2.1 Study Area and Data

Turkey is a country that is a bridge between Asia and Europe. Turkey’s geographical location is in the 36-42º North latitudes and 26-45º East longitude. The average elevation is 1130 meters above sea level and gradually increases from the West to the East of the country. The sea surrounds it on three sides (^{Mızrak, 1983}; ^{Atalay, 2011}; ^{Çelik, 2020}). Turkey has very different climatic and microclimate zones due to its geographical location and structure. Turkey is divided geographically into seven regions (Marmara, Aegean, Mediterranean, Central Anatolia, Black Sea, Eastern, and Southeastern Anatolia). It has a semi-arid climate characteristic.

Turkey has a high solar energy potential as it is close to the equator rather than the poles. According to Turkey’s Solar Energy Potential Atlas (GEPA), Turkey’s annual total sunshine duration is 2741.07 hours year^-1, its daily total sunshine duration is 7.5 hours day^-1, its daily total solar energy is 4.18 kWh m^-2 day^-1 and its annual total solar energy is 1527.46 kWh m^-2 year^-1. The Black Sea Region has the lowest radiation due to the higher number of rainy-day counts and geographical position. While Marmara and Aegean Region have radiation in medium value, Central Anatolia, East Anatolia, Mediterranean, and Southeastern Anatolia regions have high radiation levels. This study used the meteorological stations with the highest radiation (Mugla, Antalya, Mersin) and the lowest (Ordu, Trabzon, and Rize). The meteorological stations’ locations on the map of Turkey are given in Figure 1. Some characteristic features of the meteorological stations are shown in Table I.

Fig. 1 Geographical positions of stations in the study area 

Table I Some characteristic features belong to meteorological stations 

Stations	Observation Period	Elevation	Latitude	Longitude
Mugla	1969-2002	625	37.03	27.43
Antalya	1983-2006	39	37.04	31.79
Mersin	1984-2010	6	36.07	32.83
Ordu	1987-2010	4	40.98	37.87
Trabzon	1972-2005	10	41.00	39.71
Rize	1987-2019	6	41.02	40.51

2.2 Methodology

2.2.1 Empirical Equations

Several empirical equations are used in literature to explain the relationships between solar radiation and meteorological data. FAO-56 (^{Allen et al., 1998}) suggested using Angstrom and Hargreaves-Samani equations when solar radiation data cannot be measured. In the Angstrom equation, solar radiation is calculated by being associated with extraterrestrial radiation and relative sunshine duration. The Angstrom equation is given in Equation 1.

Rs= as+bs*nNRa (1)

where R _s is solar radiation at the surface (MJ m^-2 day^-1), R _a is solar radiation at the top of the atmosphere (MJ m^-2 day^-1), n/N is the relative sunshine duration (hour), and a _s and b _s are regression constant (a _s = 0,25 and b _s = 0,50).

^{Hargreaves and Samani (1982)} developed Equation 2 for solar radiation calculation by relating T_max, T_min, and R_a parameters.

Rs= kRSTmax-Tmin Ra (2)

where T _max is maximum air temperature (ºC), T _min is minimum air temperature (ºC), k _RS is adjustment coefficient (0.16- 0.19).

2.2.2 Support Vector Machine (SVM)

^{Vapnik (1995)} developed the SVM method to solve pattern recognition and classification problems (^{Cortes and Vapnik, 1995}). SVM is based on statistical learning theory and is mainly used to distinguish two data classes in the best possible way. For this purpose, decision boundaries or hyperplanes are determined. In a nonlinear dataset, SVMs cannot plot a linear hyperplane. Therefore, kernel numerals are used. The Kernel method widely increases machine learning on nonlinear data. The SVM estimator (y) process in the Kernel method is expressed as in Equation 3.

y= Kxi*Wjk+b (3)

where the Kernel core function is K _x,i , b is the bias term of the SVM network, and W _jk is called the weight vector. K _x and W represent Lagrange multipliers. K _x,i is a nonlinear function that maps input vectors to a high-dimensional feature space. Various function solutions determine the kernel core (K _x,i ). In this study, the kernel function for the SVM method was determined as polynomial due to the lack of linearity between the input and output data. In the SVM method, the ‘Batch size’ value is 100, and the ‘Complexity’ value is 0.5.

2.2.3 Artificial Neural Networks (ANN)

ANN is a mathematical modeling method whose development was inspired by biological neural systems like the human brain. In general, ANN is a system consisting of non-linear artificial cells that can be arranged as single-layer or multi-layer and work in parallel (^{Liu et al., 2010}). The basic unit of ANN is the artificial neurons. Neurons are structures that produce a response through input nodes and process information within themselves. Artificial neurons have five components: inputs, weights, combination function, activation function, and outputs (^{Cıtakoglu and Ozeren, 2021}). Artificial neurons get together to create ANN. The most common structure of ANN consists of the input layer, where the data is presented to the artificial neural network as input data; the hidden layer, where the data is processed; and the output layer, where the result is obtained. The basic structure of ANN is given in Figure 2.

Fig. 2 Artificial Neural Networks Architecture 

ANN has been the most important tool in solving complex nonlinear problems. ANN is widely used in hydrological studies such as precipitation-flow models, precipitation, runoff, temperature, snowmelt, evapotranspiration, solar radiation, and drought. (^{Khotanzad et al, 1996}; ^{Dawson and Wilby, 1998}; ^{Machado et al, 2011}). This study used four hidden layers for the ANN model, 20 neurons in each hidden layer, the hyperbolic tangent function as the activation function, and the Levenberg-Marquardt algorithm as the training algorithm.

2.2.4 Long Short-Term Memory Networks (LSTM)

LSTM, a type of recurrent neural network (RNN) approach, was first proposed by ^{Hochreiter and Schmidhuber (1997a)}. Many people then developed the method. The LSTM model was specifically designed to overcome the exploding and disappearing gradient problems that typically arise when learning long-term dependencies, even if the minimum time delays are very long (^{Hochreiter and Schmidhuber, 1997b}; ^{Van Houdt et al., 2020}). The LSTM unit consists of a cell, an input gate, an exit gate, and a forget gate. The structure of LSTM is given in Figure 3. These gates handle the writing, reading, and resetting of the cell. The cell remembers values at random time intervals, and three gates regulate the information flow associated with the cell (^{Singh et al., 2017}). In the LSTM network structure, the input data is multiplied by the output of the input gate to identify new information that can be accumulated in the cell. To calculate the information that can be propagated to the network, the output data of the network is multiplied by the activation of the output gate. The cell states of the previous time are multiplied by the activation of the forget gate to determine whether the last state of the cell is forgotten. (^{Sainath et al., 2015}).

Fig. 3 LSTM Block (^{Sagheer and Kotb, 2019}) 

(X _t ) is an input at time step t, the hidden state from the previous time step (S _t-1 ) is introduced to the LSTM block, and the hidden state (S _t ) is then calculated as follows:

The first step; is to decide what information to discard from the cell state. This decision is made by the following forget gate (ft):

ft= σ( Xt*Uf+ St-1*Wf+bf) (4)

The second step; is to decide what new information to store in the cell state. There are two stages in this step.

1. The input gate (i _t ) layer decides which values to update.

2. The tanh function determines the candidate information that will create new information.

it= σ(Xt*Ui+ St-1*Wi+bi) (5)

Ćt=tanh(Xt*Uc+ St-1*Wc+bc)  (6)

The third step, from the old cell state C _t−1 is to create the new cell state Ct as follows.

Ct= Ct-1*ft+ it*Ćt) (7)

The final step decides what will be produced as output. The output gate (O _t ) determines which components of the cell state will be generated as output in this stage. Then the cell state goes through the tanh layer (to push between -1 and 1 the values) and multiplies it with the output gate as follows

Ot= σ(Xt*U0+ St-1*W0+bo)  (8)

St= Ot*tanh⁡(Ct) (9)

From the previous six equations, the LSTM presents the following three groups of parameters:

1. Input weights: U ^f , U ⁱ , U ^c , U ⁰ .

2. Recurrent weights: W ^f , W ⁱ , W ^c , W ⁰

3. Bias: b _f , b _c , b _i , b _o .

The LSTM model used in the study has one input layer, one hidden layer containing 64 LSTM neurons, and one output layer that makes a single value estimation. Sigmoid was used as the activation function in LSTM neurons. The learning rate of the model was determined as 0.005. The dropout value is set to 0.1 to avoid over-fitting. The stochastic Gradient Descent (SGD) algorithm was used as the optimization method. The LSTM model is trained for 100 epochs, and one batch size is used in this training (^{Li and Cao, 2018}; ^{Granata and Di Nunno, 2021}).

3. Results and Discussion

The average MAE, RMSE, MRE, and R ² values obtained from the experiments performed for the Hargreaves-Samani and the Angstrom equation are given in Tables II and III.

As can be seen in Tables II and III, two different models created with input from the Hargreaves-Samani and Angstrom equations for the LSTM model are the most successful among machine learning and deep learning methods used in experimental studies. The results obtained with machine learning and deep learning models show more successful outcomes than the solar radiation values calculated by both the Hargreaves-Samani and the Angstrom equations alone. The results obtained with input from the Angstrom equation are more successful than the results with input from the Hargreaves-Samani equation. MAE, RMSE, and MRE values approaching zero and R ² values approaching one indicate the success of our proposed model. The graphical representation of the absolute errors in the test dataset for solar radiation value estimation using the LSTM deep learning model, for the most successful model with the Hargreaves-Samani equation, is shown in Figure 4.

Fig. 4 Graphical representation of absolute errors of the Hargreaves-Samani equation and LSTM model 

As seen in Figure 4, the Hargreaves-Samani equation and LSTM deep learning model were compared according to the monthly solar radiation value estimation results with a dataset consisting of 2112 records and eight attributes. The LSTM deep learning model was observed to be more successful than the Hargreaves-Samani equation and minimizes the estimation error. It was observed that the MAE value, which was 2.65 with the Hargreaves-Samani equation, decreased to 0.987 with the LSTM model.

The graphical representation of the absolute errors in the test dataset for solar radiation value estimation using the LSTM deep learning model, for the most successful model with the Angstrom equation, is shown in Figure 5.

Fig. 5 Graphical representation of absolute errors of the Angstrom equation and LSTM model 

As seen in Figure 5, the LSTM deep learning model was also more successful than the Angstrom equation and minimized the estimation error. As a result of the experimental studies, although linearity is observed between the solar radiation values and the input data set, it was observed that the error value was reduced by using SVM and ANN. In this situation, the LSTM model, one of the deep learning methods, was applied, and it was observed that the MAE value, which was 1.24, decreased to 0.747 with the Angstrom equation. In addition, although the solar radiation estimation calculated using the Angstrom equation seemed more successful than the Hargreaves-Samani equation, the LSTM model showed more successful results than both equations. The test samples shown in Figures 4 and 5 were randomly taken from six different stations. The random sampling of the test samples showed that the data collection periods were also random. This is especially important in terms of objective evaluation of the results.

Figure 6 and Figure 7 showed box-whisker error plots corresponding to the absolute values of the difference between solar radiation values calculated from the Hargreaves-Samani and Angstrom equations and estimated values from the models. In this way, the distribution characteristics of the errors were revealed, and the models were compared objectively (^{Di Nunno et al., 2021}). The line in the middle of the box in the box-whisker plot shows the median value (Q ₂ ), the lower boundary line of the box indicates the first quartile (Q ₁ ), and the upper boundary line of the box indicates the third quartile (Q ₃ ). As the Q ₃ - Q ₁ (IQR-interquartile range) value increases, the spread of the data increases, and as it gets smaller, it decreases. As the Q ₂ - Q ₁ value gets larger, the median value moves away from the first quarter, and the data set shows a left-skewed distribution. As the Q ₃ - Q ₂ value gets larger, the median value moves away from the third quartile, and the data set shows a right-skewed distribution. The data set exhibits a symmetrical distribution when the median, mean, and mode values are in the middle of the quadrants. The cross on the box gives the mean of the data set. Therefore, it shows the MAE value. The vertical lines coming out of the boxes denote whiskers. The lower whisker or lower line represents the smallest (minimum), and the upper whisker or upper line represents the largest (maximum) data value. The dot or dots outside the whiskers represent outliers. The position of the quadrants changes as the number of outliers and their spread beyond the whiskers increases. Therefore, when comparing multiple box-whisker plots, it is necessary to consider the number of outliers, their spread, the box width, and the whiskers’ lengths.

Fig. 6 Boxplot of Hargreaves-Samani, ANN, LSTM, SVM absolute errors 

Fig. 7 Boxplot of Angstrom, ANN, LSTM, SVM absolute errors 

Figure 6 illustrates that while the Hargreaves-Samani model does not contain outliers, the ANN, LSTM, and SVM models do. Although the Hargreaves-Samani model appears to be more favorable in this regard, it is clear that the LSTM model stands out when compared to the others when looking at whisker lengths, box sizes, MAE values, and how far these values are from the quarters. On the other hand, the Hargreaves-Samani model exhibits a right-skewed structure, and the SVM model exhibits a left-skewed structure. ANN and LSTM models show a more symmetrical structure. Considering the whisker lengths, LSTM can be preferred over the ANN model.

As seen in Figure 7, although all models show a right-skewed distribution, the size of the boxes and whisker lengths in the Angstrom model showed a wider spread than the ANN, LSTM, and SVM models. Although the box widths of the ANN, LSTM, and SVM models seem the same, the smallest data value in the LSTM model almost coincided with the Q1 quartile. This situation caused the MAE value to be even smaller. Therefore, in terms of the distribution of absolute error values, the most stable model seems to be LSTM.

4. Conclusion

This study has developed LSTM deep learning, ANN, and SVM machine learning models to estimate Turkey’s monthly solar radiation value. Training and testing processes were carried out with the meteorological data of a total of six stations in Turkey that receive the most and the least solar radiation. When the developed and empirical models are evaluated together, it is concluded that all learning models make more successful estimations than empirical models. While the LSTM model has lower MAE, RMSE, and MRE values than other learning models and empirical models, it has a higher R ² value. Therefore, the best-performing method is the LSTM deep learning model. As a result, experts who perform solar radiation calculations can successfully estimate solar radiation with the LSTM model in cases where measurement is impossible.

This study, which deals with estimating the solar radiation value, transfers the computational engineering task from the researchers to the deep learning model. The deep learning approach performs much better when comparing the deep learning method, the empirical models, and the traditional machine learning models. The success of the deep learning approach is based on two main reasons. The first is that deep learning models mimic how the human brain works. High-level attributes that represent semantic relationships are extracted by taking the lower-level attributes with the successive layers of the model. Second, thanks to the memory module in the LSTM model, it automatically identifies active features at each step of the training. Finally, for the estimation process used in this study, ANN, SVM, and LSTM machine learning/deep learning methods seem to be methods that can be used in relational screening studies. These methods are known to be powerful estimation methods. For this reason, it is recommended that these methods be used more widely in similar studies to be conducted in the future.