Parameter Optimization of Rainfall-runoff Model GR4J using Particle Swarm Optimization on Planting Calendar

Challenges in food production in the future are certainly more complex in developing countries like Indonesia. The Agricultural Research and Development Agency developed an Integrated Planting Calendar (Katam) information system to decide the water discharge for rice planting time using rainfall-runoff modeling called GR4J (Genie Rural a four parameters Journalier). This study aims to improve the accuracy of the GR4J model for determining water discharge in an area. The study areas in this study were the Progo Watershed and the Wuryantoro Watershed. The GR4J model is measured based on four free parameters in the form of Maximum Capacity Production Store (X1), Coefficient Changes in Groundwater (X2), Maximum Capacity Routing Store (X3), and Peak Time Ordinate unit hydrograph (X4). The four parameters are optimized using Particle Swarm Optimization (PSO). This study shows that the parameter optimization of the GR4J model with PSO was successfully carried out. As a determinant of the model's success, the Nash-Sutcliffe Efficiency (NSE) equation and the Root Mean Square Error (RMSE) is used. We have found that NSE is increasing and RMSE decreases in RMSE in each watershed after using PSO optimization. The area affects the accuracy of the GR4J model. The smaller the area, the more it can show the characteristics of a watershed. Using the PSO algorithm in the GR4J model will be more effective than using manual or trial methods based on expert knowledge because it can achieve optimization quickly. Keywords— GR4J model; optimization; particle swarm optimization; planting calendar.


I. INTRODUCTION
Rice is one of the significant agrarian items in Indonesia [1]. More than half of Indonesia's population depends on rice produced from rice plants [2]. Rice is defenseless against environmental change with the goal that environmental change influences rice creation and even neglects to reap. Because of the dry season and overflowed, many harmed rice manors. Assurance of the growing season can answer environmental change transformation, so rice creation stays stable [3]. In recent decades, food crop farming in Indonesia often relies solely on hereditary habits in determining the planting time. As a result, farmers are often faced with water shortages, especially in high rainfall intensity over a short period or an extended dry period [4]. Local wisdom and conventional methods used to apply cropping patterns have experienced a bias due to early shifts in the growing season. Determination of the right planting time is one of the determinants of crop success and the recall of food crop productivity. One way to overcome water limitations by using rainfall data, river water, and irrigation water needs is to determine planting time [5]. The Schedule of planting individual plants for a year in an area, including land preparation, planting, and harvest, is called the planting calendar [6]. Indonesia has an information system for determining planting rice, namely the Integrated Planting Calendar (KATAM). KATAM provides information on estimating the beginning of planting time, the potential for planting area, crop rotation, and cropping intensity at the sub-district level for one year [7]. KATAM was developed by the Research Institute for Agro-climate and Hydrology Indonesia (Balitklimat) in determining the water discharge for rice planting time using rainfall-runoff modeling called GR4J (Genie Rural a four parameters Journalier).
The GR4J model is estimated dependent on four free parameters; they are the maximum capacity of the production store (X 1 ), groundwater exchange coefficient (X 2 ), maximum capacity of routing store (X 3 ), and time peak ordinate of hydrograph unit UH1 (X 4 ) [8]. The GR4J model was developed from the previous model, namely GR3J (Genie Rural three parameters Journalier), which was developed from the original method found by Edijatno and Michel, then it was refined by Nascimento and Edijatno [9]- [11]. This model has subsequently been applied to various watersheds such as Qu River Basin in East China, five basins in northern Tunisia, the Lake Tana Basin in Ethiopia, Heshui watershed in Jiangxi province, and Songhua River Basin in northeastern China [12]- [16]. The optimization of the GR4J model parameters has been carried out using genetic and Gauss-Newton algorithms. Based on these studies' results, genetic algorithms and Gauss-Newton successfully increased Nash criteria with values above 70% from the previous value of only 48% [17].
Although the GR4J model parameters' optimization has been done with genetic algorithms and the Gauss-Newton method, research that focuses on optimizing the GR4J model parameters with other optimization algorithms, such as Particle Swarm Optimization (PSO), has never been done. Among the many techniques that can be used as optimization tools, PSO is one of the newest and efficient optimization tools [18]. PSO calculation was planned as a straightforward, simple, and practical calculation to tackle a few streamlining issues [19]. Whereas in some studies, PSO algorithms generally have performance that is as good as genetic algorithms. Even in some cases, PSO can exceed the performance of genetic algorithms [20]. Another advantage is that PSO has fewer function evaluations than genetic algorithms, so that it has a relatively better level of computational efficiency. Based on these considerations, this study examines aims to optimize the parameter value of the GR4J model with the PSO algorithm.
This study examines the Progo watershed and Wuryantoro sub-watershed to provide an overview of the simulated water discharge that will be carried out. Based on previous studies' results, in this experiment, we will use the PSO algorithm to obtain the optimal value of each free parameter of the GR4J model in determining the simulation water discharge. The results of the model will be used as a recommendation for determining KATAM.

II. MATERIAL AND METHOD
The study areas are the Progo Watershed and the Wuryantoro sub-Watershed in Indonesia. The Progo watershed has an area of around 1762.7 km 2, and The Wuryantoro sub-watershed has an area of about 6 km 2 . The Progo watershed is a watershed area across Central Java Province and Special Region of Yogyakarta, Indonesia, and the Wuryantoro sub-watershed located in Wonogiri Regency, Central Java, Indonesia.
The data used in the Progo watershed and the Wuryantoro sub-watershed consist of rainfall and daily observation discharge in 2006, which were collected from the Research Institute for Agro-climate and Hydrology. The data used to calculate the potential evapotranspiration of the Progo watershed and the Wuryantoro sub-watershed include geographical location, air temperature, wind speed, air humidity, and solar radiation, and altitude. The data were daily data in 2006 and can be downloaded via the site http://ogimet.com/gsynres.phtml.en. Ogimet is a Weather Information Service. Ogimet was created by Guillermo Ballester Valor, a meteorologist from the Spanish Meteorological Institute, in 2005. His contribution was also obtained from Dr D. José Manuel García López, Catedratic from the Department of Histology of Granada University. Distributed information is utilized by a few world data systems, particularly from the National Oceanic and Atmospheric Administration (NOAA), and uses free permit programming to deal with it. The Research Institute currently uses Ogimet data for Agroclimate and Hydrology to calculate potential evapotranspiration. The research stages are illustrated in Fig. 1.

A. Determination of Potential Evapotranspiration Values
Potential evapotranspiration (E) is the process of the occurrence of maximum water loss from a surface covered by short grass to the atmosphere, which is influenced by meteorological factors [21]. Daily potential evapotranspiration analysis in 2006 was calculated using Microsoft Excel software. In calculating potential evapotranspiration, climatology data are needed, including geographical location, air temperature, wind speed, air humidity, solar radiation, and altitude (altitude). Microsoft Excel software for calculating potential evapotranspiration in this study comes from the Indonesian Institute of Sciences calculations. The calculation of potential evapotranspiration uses the Penman Modification method, measured based on Equation 1 [22].
where ETP represents evapotranspiration (mm/day). W represents the weighting factor. R n represents net radiation, equivalent to evaporation (mm/day). f (U) represents the function of wind speed. e s -e a represents the contrast between saturated vapor pressure and actual vapor pressure on average at the average temperature (mbar). c is the adjustment factor to compensate for the influence of weather conditions during the day and night.

B. GR4J Model Process
The daily debit model applied in the Integrated Katam information system is GR4J. This model is a simple, reliable debit model and is a global model with a daily time interval that only requires four parameters [23]. This research processes the GR4J model using the package airGR in R software [24], [25]. Daily rainfall (P) and potential evapotranspiration (E) are inputs for running the GR4J model. The observation flow discharge data (QMi) are used to observe the comparison of the results of the simulation flow discharge (QSi) with the observation flow discharge (QMi) used as a determinant of the success of the model. The GR4J model has four parameters, namely the maximum capacity of production deposits (X 1 ) in mm units, groundwater exchange coefficient (X 2 ) in mm units, maximum capacity of the routing store (X 3 ) in mm units, and the peak time of the hydrograph unit (UH) (X 4 ) in units of days. Figure 2 indicates the process of the GR4J model. Figure 2 explains that the GR4J model uses the input of daily rainfall (P) and potential evaporation (E). Rainfall falls on the surface and fills soil moisture; the volume is indicated by parameter X 1 . However, direct runoff can occur if the intensity of rain is high. Deposits at X 1 undergo a flow process and meet with direct runoff, and finally, 90% of the flow propagates to deeper deposits where the volume is indicated by parameter X 3 . The remaining 10% from the previous process produces runoff flowing to the outlet set with a unit hydrograph in the parameter. Deposits inside are also still propagating and produce delayed runoff. Both river deposits and runoff are susceptible to internal infiltration outside the system in exchange for parameter X 2 . The output of the GR4J model is the simulated flow discharge (QSi). Parameters X 1 , X 2 , X 3 , and X 4 , are parameters that affect the simulation flow debit (QSi). Determination of parameters X 1 , X 2 , X 3 , and X 4 by the Indonesian Institute of Sciences is conducted randomly until the simulation flow rate (QSi) is obtained, which approaches the observation flow's discharge flow (QMi). In this study, the experimental process for determining parameters X 1 , X 2 , X 3 , and X 4 was carried out using the PSO algorithm.

C. GR4J Parameter Optimization with PSO algorithm
Particle Swarm Optimization (PSO) is a ground-breaking nature-roused metaheuristic improvement strategy. Contrasted with different strategies, PSO can decide the ideal arrangement in fewer assessments and, for the most part, performs even more efficiently and effectively [26]. The initialization of the PSO algorithm system usually starts in a random population. The population is usually called particles. Each particle is a potential solution in an optimization problem. Particles are still heading towards the optimal solution possible. The speed of movement is influenced by the velocity at which each iteration is repeated. Changes in each particle's velocity are determined by previous velocity, pbest location, and gbest values. Various random values are generated as acceleration pbest and gbest [27]. Pbest (personal best) is a particle that gives a better solution than the previous particle, while gbest (global best) is the best solution that particles have ever obtained in a population. Fig. 2 The GR4J rain-flow surface model [28] PSO is a population-based search algorithm. PSO is an algorithm based on the swarming behavior of animals in nature [29]. In PSO, particles can travel randomly in any direction possible towards their optimal global and actual best position using their present velocity, modified later than randomly. Thus, the use of real-number randomness and global particle communication can increase swarm efficiency [30]. In the parameter optimization stage, using the package from hydroPSO in software R. Optimization is done on the parameter values of the GR4J model, namely X 1 , X 2 , X 3 , and X 4 . The results of parameter optimization with the PSO algorithm are used to determine the optimal parameter values based on the evaluation of the GR4J model. The stages of the PSO method in this study are as follows: 1. Initialize each parameter of the GR4J model, namely X 1 , X 2 , X 3 , and X 4 as populations, and each population has particles with random position and speed in dimension d in problem space. The problem space in the GR4J model is assumed to be the search limit, which is composed of the upper and lower limits. 2. The test uses a fitness function for each particle. In this case, the fitness function is an optimized function of the Nash-Sutcliffe Efficiency (NSE) model. The suitability of the model is better if NSE approaches 1 (positively correlated). 3. Compare the particles evaluated earlier with the pbest (personal best) particles. If the value is greater than pbest, then the pbest value set is equal to that value, and the pbest position is the same as the particle location evaluated using the fitness function.
4. Identify the best particles from a population as gbest (global best). 5. Check the particle velocity and position with Equations 2 and 3 [26].
Where V id represents dimensional particle velocity, and X id represents the position of particle dimension currently. C 1 is the acceleration coefficient of personal influence, and C 2 is the acceleration coefficient of social influence. P id represents pbest, and P gd represents gbest. 6. Repeat the second step until the stop criteria are met, usually reaching the optimal value, or reaching a certain number of iterations.

D. Evaluation of the GR4J Model
The GR4J model's performance was evaluated using Nash-Sutcliffe Efficiency (NSE) and Root Mean Square Error (RMSE). NSE can be used as a benchmark for the success of the GR4J model. The NSE general equation is formulated as in Equation 4 [31].
Where , is discharge observation,!̅ / is average observation debit value, 1 is simulation model discharge, and 0 3 / is the average simulation model debit value. NSE has a range between 0 and 1. NSE has several criteria, as shown in Table 1 [32].
where RMSE is Root Mean Square Error, Q Si is simulation model value, Q Mi is observation value, and n is the amount of data.

A. Determination of Potential Evapotranspiration Values
Data used to calculate potential evapotranspiration in the Progo watershed and Wuryantoro sub-watershed consist of geographical location (latitude), temperature, humidity, wind speed, solar radiation, and day length. The data is daily data in 2006 for the Progo watershed and the Wuryantoro subwatershed. Data can be downloaded via the site http://ogimet.com/gsynres. phtml.en. Ogimet is a Weather Information Service. Ogimet was created by the Spanish Meteorological Institute meteorologist Guillermo Ballester Valor. Dr D. José Manuel García López, Catedratic from Granada University's Department of Histology, has also contributed to Ogimet 's development. Published data is also collected by many world databases networks, mostly by the National Oceanic and Atmospheric Administration (NOAA), stored using free license software. Ogimet data is data currently used by Balitklimat to calculate potential evapotranspiration. The dataset attributes used for evapotranspiration calculations and the GR4J model can be seen in Table 2.   TABLE II  The Progo watershed results and the Wuryantoro subwatershed using Equation 1 can be seen in Fig. 3. The highest potential evapotranspiration in the Progo watershed occurs in November, and the lowest occurs in May and December. The potential evapotranspiration rate ranges from 1.56 mm/day to 5.33 mm/day, and the average is 3.54 mm/day. In the Wuryantoro watershed, the potential evapotranspiration was highest in March and April, while the lowest is formed in January. Potential evapotranspiration rate ranges from 2.4 mm/day to 3.5 mm/day and an average of 2.99 mm/day. This stage aims to obtain optimal values from parameters X 1 , X 2 , X 3 , and X 4 in the GR4J model with the PSO algorithm. Searching for parameter values X 1 , X 2 , X 3 , and X 4 is done randomly, with search limits X 1 , X 2 , X 3 , and X 4 = {-9.99 ... 9.99}. After the value of the parameters X 1 , X 2 , X 3 , and X 4 was determined, the next step is to process the maximum GR4J model as many as 1000 repetitions until it gets the optimal parameter value. The exponential function gets the value of the transformation from the initial value using the optimal parameter values X 1 and X 3 . As for the parameter value, X 2 entered the SINH function to return the initial value's hyperbolic sine value. Parameter X 4 value was also entered into an exponential function. However, it was added by 0.5 to get the value of the transformation from the initial value. The initial values' results before transformation and after transformation in parameters X 1 through X 4 optimized by Table 3 displays the PSO algorithm.  The evaluation was carried out on each loop of the GR4J model by optimizing parameters using the PSO algorithm until an NSE value approaching 1 and RMSE approaching 0. The evaluation was carried out by implementing Equation 4 to calculate NSE and Equation 7 for RMSE calculation. The best evaluation results for the Progo watershed and Wuryantoro sub-watershed is given in Table 4.

C. Evaluation of the GR4J Model
The GR4J model process results with optimal parameters for simulation debit modeling Progo watershed and Wuryantoro watershed can be seen in Fig. 4. In Fig. 4, there is a graph with a black line, which is a line to indicate the observational discharge taken from the discharge data in each location, while the red line is a simulation discharge originating from the GR4J and blue line modeling which is a daily rain derived from the daily rainfall data of each watershed.
In Fig. 4, the flow rate of the simulation results with the observation flow discharge has a similar pattern. Based on The maximum water discharge occurs on January 5 at 11.64 mm/hour, and on October 22 a minimum water discharge of 0.28 mm/hour occurs. Based on Fig. 4 (b), the daily average daily discharge of the Wuryantoro subwatershed in 2006 is 2.36 mm/hour. The maximum water discharge occurs on December 28 at 58.59 mm/hour, and during September, a minimum discharge of 0.14 mm/hour occurs.
Based on the results of the study, it showed the relationship of simulated discharge (Qsi) with discharge observation (QMi). In Fig. 5, it can be seen that the scatter plot is adjacent to the slope line (slope = 1), which explains that the simulation discharge has a strong relationship to the observation discharge. The NSE is 0.7188 in the Progo watershed, which is categorized as strong, and the NSE is 0.8632 in the Wuryantoro sub-watershed, which is categorized as very strong. In addition to calculating the NSE, an error calculation from the model with RMSE was also obtained. The RMSE is 1.3944 in the Progo River Basin, and the RMSE is 2.0055 in the Wuryantoro sub-watershed.

D. Determination of KATAM
Determination of KATAM in an area is done by the analysis of rainfall and water discharge data. After getting the optimal parameters in each region, these parameters will be used to determine water flow in the following year. Balitklimat collaborates with the Meteorology, Climatology, and Geophysics Agency (BMKG) to obtain rainfall forecast data. Rainfall forecast data is processed into the GR4J model by entering optimal parameters in order to get the water discharge forecast data. Determination of KATAM in an area is not only based on the approach using the GR4J model. Other approaches such as considering natural and social conditions as well as resource-based, for example, using agricultural statistical data, using the assumptions and subjectivity of experts. Analysis of the GR4J model used to provide a general description of water availability in the Progo watershed and Wuryantoro sub-watershed has been successfully carried out. Optimizing the GR4J model parameter values with the PSO algorithm has produced the best evaluation values with optimal parameters. This value can be seen from the increase in NSE and decreased RMSE in each watershed after using PSO optimization. The area size affects the accuracy of the GR4J model. The smaller the area, the more it can show the characteristics of a watershed. Progo watershed produces an NSE value of 0.7188 and RMSE of 1.3944, indicating that the Progo watershed is categorized as vital. The Wuryantoro sub-watershed produces an NSE of 0.8632 and RMSE of 2.0055, indicating that the Wuryantoro sub-watershed is categorized as very strong. The PSO algorithm used in the GR4J model will be more effective than using manual methods or trial and error based on expert knowledge because it can achieve optimization in a relatively short amount of time. The Agricultural Research and Development Agency can optimize GR4J model parameters with PSO to develop a KATAM information system.