Multi-objective optimization of groundwater monitoring network using a probability Pareto genetic algorithm and entropy method (case study: Silakhor plain)

Hydrometric monitoring networks play an important role in water resources planning and management. Groundwater monitoring network design can be divided into two streams: the stream designed for groundwater quality monitoring and that designed for groundwater quantity monitoring. This does not mean that one network could not perform both tasks, but it means that the approaches used in the design vary based on the intended use of the network. Hence, a suitable monitoring network needs to maintain and expand due to the importance of hydrometric data for water resources management practices as the availability and reliability of hydrometric data are important for planning and management purposes (Mishra & Coulibaly 2009).

Groundwater monitoring approaches have two main categories, hydrogeological and statistical simulation, variance or probability-based methods (Loaiciga et al. 1992). For optimal groundwater management, it is necessary to collect sufficient information from the quantitative and qualitative features of the aquifer. Assessing and anticipating the level of groundwater through specific models help to predict groundwater resources. The design of observation networks is a multi-objective optimization problem. In the case of multi-objective optimization, measurements are highly considered in terms of spatial and temporal frequencies. On the other hand, these calculations require a large amount of infrequent data intervals, resulting in costs going up significantly. To this end, the network must be optimized to a degree that adequately reflects the level of the station in the target area. In order to generalize the measurement results from available wells to other points without statistics, the well network should have a density appropriate to the area (Hosseini & Kerachian 2017; Keum et al. 2017; Ayvaz & Elçi 2018). To implement this, different models and algorithms such as land-based methods and optimization algorithms have been used. In recent years, with the demonstration of the capabilities of intelligent models, these models have been applied in groundwater modeling. Intelligent models such as artificial neural networks (ANN) and support vector machines (SVM) are more developed than other methods due to their processing speed, high accuracy, and relatively low cost. There have been a wide range of studies on the design of groundwater monitoring networks (e.g., Bashi-Azghadi & Kerachian 2010; Kollat et al. 2011; Barca et al. 2013; Reed & Kollat 2013; Khader & McKee 2014; Bhat et al. 2015; Júnez-Ferreira et al. 2016; Komasi et al. 2018; Pourshahabi et al. 2018). In this way, Asefa et al. (2004) used the SVM method based on Statistical Learning Theory (SLT) to design a long-term groundwater monitoring network to reduce spatial redundancy. In this study, the possibility of generalizing the developed method for the analysis of time redundancy has been mentioned. Another approach in terms of groundwater monitoring networks is entropy theory. In this way, Mogheir et al. (2003) used entropy theory in research to design a site for sampling groundwater quality monitoring networks. They defined the collection of groundwater quality variables using the information transfer model and the correlation model which determine the location transport model of the wells better than the correlation model. In another study, Leach et al. (2016) considered feeding groundwater as the basis for the design of a monitoring network using the dual entropy optimization model to design an optimal groundwater monitoring network. Utilizing entropy is a suitable method for finding common information among observation wells for optimization. Some studies have a scattered method for designing a groundwater monitoring network. Li & Hilton (2005) used the method of optimizing ants to solve the problem of optimizing the design of the groundwater monitoring network. In this study, while investigating a pollutant in the quality monitoring network, only a spatial optimization has been achieved. Also, Guo et al. (2011) compared several methods, sampling, hydrological analysis, and a method of non-homogeneous mean level approach to optimizing the observation network of groundwater in the field of designing a network for observing groundwater level, which results in the superiority of this method over other methods. In another study, Bhat et al. (2012) used the usual Kriging method to design a groundwater monitoring network in Florida. The results showed that 58 wells should be reduced to 44 wells in the region to reduce the estimation error. Zhou et al. (2013) utilized mapping to optimize the observation network for groundwater level. The results show that these maps can distinguish the regions with unique spatial and temporal characteristics, and can be used to design the grid. Also, Esquivel et al. (2015) used multi-criteria analysis to design an optimal groundwater table monitoring network. A multi-criteria analysis has been implemented with the help of GIS on an aquifer in Mexico. The analytic hierarchy process (AHP) method was used to weigh the criteria and the stability rate was 0.08.

Multi-objective optimization algorithms have been used in network design. For example, Kollat & Reed (2007) used the epsilon-dominance non-dominated sorted genetic algorithm II (NSGA-II) to design long-term groundwater monitoring networks based on cost, contaminant estimation error, uncertainty, and contaminant mass error. Later, Kollat et al. (2008) studied the next generation epsilon-dominance hierarchical Bayesian optimization algorithm (ε-hBOA) and its application to groundwater network design and found that it was superior for network design due to its model building capability. In a further study, Dhar & Patil (2012) presented a method for designing a groundwater quality monitoring network under cognitive uncertainty conditions in predicting the mass flow through the second version of the genetic algorithm (GA) and Kriging. The results of this study showed the effectiveness of the proposed method with two approaches. Leach et al. (2015) added objective functions for stream flow signatures and indicators of hydrologic alteration to explicitly consider the shape of the hydrograph at different gauging sites when building a hydrometric network. This modification increases the variety of flow regimes seen in optimal networks. Luo et al. (2016) aimed at optimizing the long-term groundwater monitoring networks using PPGA. In this research, four objective functions were considered for optimization of the observation network. Using the PPGA, efforts to obtain the Pareto-front of four targets have been considered.

According to the literature review, the present study aims to assess the optimal groundwater monitoring network with the application of GA and Kriging as interpolation methods regarding the case study while focusing on the Silakhor plain located in Iran. This study designs an optimal groundwater network in two scenarios. The first scenario redesigns and evaluates the groundwater monitoring network by entropy theory, and the second scenario chooses optimal wells from the existing wells by considering entropy of each station in the study area concerning the groundwater monitoring network.

The data used in this paper are from the Silakhor plain located in the west of Iran in Lorestan province. The Silakhor plain is an upper district of Karun River basin. Starting from Garrin Mountain, Gelerood River passes southern Borujerd and receives more creeks in lower Silakhor where it is called the Tireh River and is an artesian aquifer. The Silakhor plain is located between 33° 15′ 46″ and 34° 10′ 00″ northern latitude and 48° 28′ 53″ and 49° 30′ 25″ eastern longitude. The Silakhor aquifer area is 2,545.8 km² (Figure 1). For modeling and designing the optimal monitoring network, the statistical data of the Silakhor plain for a ten-year period were used. The statistical characteristics of its stations are presented in Table 1.

Figure 2 shows monthly rainfall time series and groundwater level changes in the Silakhor plain from 2003 to 2013. There are 29 different stations in this plain. The aquifer's parameter is measured by a piezometric well in each area.

Kriging is a class of statistical techniques for optimal spatial prediction. It has been used in many other disciplines, including agriculture, mining, and environmental sciences. Kriging is a probabilistic predictor assumed as a statistical model for the data. Kriging predictors have standard errors that quantify the uncertainty associated with the predicted values. Kriging predictors are called optimal predictors because the prediction error is minimized and, on average, the predicted value and the true one coincide. Kriging uses a semivariogram (a function of the distance and direction separating two locations) to quantify the spatial dependence in data (Krivoruchko 2012).

The IDW is one of the most applied and deterministic interpolation techniques in the field of hydrology. IDW estimations were made based on nearby known locations. The weights assigned to the interpolating points are the inverse of their distance from the interpolation point. Consequently, the close points are made up to have more weights (more impact) than distant points and vice versa. The known sample points are implicit to be self-governing from each other (Preziosi et al. 2013).

Genetic algorithms can find optimum results in a multi-objective optimization problem. In this way in the GA, unlike other classical methods, a random population solution is selected. Each solution for the problem is represented as a set of parameters which are known as genes. Joining genes creates a binary bit string of values, denoting each member of the population referred to as a chromosome. A chromosome evolves through iterations, called generation (Gen et al. 1997).

To solve a multi-objective problem, two different strategies can be applied either by aggregating the weighted objectives to form a single objective problem or by finding the multiple solutions on a Pareto front to generate the best alternatives. The first method provides the leverage to solve the problem as a single objective, but assigning weights can be challenging for most problems. The second method requires solving the problem for all the objectives to obtain the non-dominated or Pareto set of solutions. Non-dominated sets of solutions consist of feasible optimal solutions (Puri et al. 2017).

In this study, monthly observation groundwater levels of the Silakhor plain from 2003 to 2013 were used in the optimal model. Using the Kriging method and time series data obtained from observation wells, potential wells were interpolated in ten years by Arc GIS. This information is considered as observation data with simulation data. Simulation data have been produced by the IDW method for ten years in potential wells by MATLAB software. According to Equation (4), the second step for determining the optimal groundwater monitoring network is to calculate RMSE between observation and simulation values for each well as the objective function. After these steps, the optimization algorithm starts to solve the optimization problem and find the optimal monitoring network. The optimization process is illustrated in Figure 3. Results of interpolation values in potential wells during ten years are shown in Figure 4.

In this research, two scenarios were examined: (1) redesigning the groundwater monitoring network in the Silakhor plain and (2) choosing the optimal wells from the existing ones in the study area for the groundwater monitoring network. Multi-objective optimization design problems of the groundwater monitoring network were solved by NSGA-II. The NSGA-II parameters include number of iterations, total population, crossover percent, and mutation percent given as 200, 1,000, 0.7, and 0.2, respectively.

In the first scenario, 29 wells were used and interpolated as observation wells in the Silakhor plain to determine the groundwater level on potential wells. Figure 5 shows the existing and proposed potential wells in the Silakhor plain. To specify potential wells’ locations, ASTM standard D 5092-90 (1996) was used. In this regulation, a grid of 492 points with 1 km distance between each point was used. A groundwater monitoring network was designed with potential wells. After determining groundwater level in potential wells as observation data, the simulated groundwater level was calculated by the IDW method in Matlab software. A matrix was introduced into the optimization model of 492 rows and 120 columns presenting the number of potential stations and monthly groundwater level in the potential stations, respectively. Then, NSGA-II looks for the optimal stations among the 492 potential stations and specifies the number of optimal stations.

According to the results of the calculations conducted by NAGA-II in Matlab software, the optimal monitoring networks were obtained with 12 observation stations. Since the response obtained from the optimizer algorithm consists of a combination of the number of stations and their position, it maps the response optimizer algorithm in the form of a Pareto-front which represents the number of optimal stations versus RMSE. The Pareto-front for the first scenario is shown in Figure 6.

As shown in Figure 6, 12 optimal wells have RMSE value of 0.61 m. Figure 7 shows the spatial distribution of 12 optimal monitoring stations. As shown in Figure 7, the program identifies a monitoring network with 12 stations, which has reduced the number of stations by 60% as compared to the existing monitoring network (i.e., 29 observation wells). The RMSE corresponding to the number of 29 stations is 0.54 m. At first glance, it can be concluded that increasing the number of monitoring stations will increase network accuracy and reduce the value of RMSE. But the other objective function for minimizing network monitoring costs avoids network recommendation with more stations because increasing the number of stations will increase the costs of building and monitoring the network. Comparing the network with 12 stations with the existing network (29 stations) shows more stations (17 stations) could reduce only 12% of the network error. Another important point is that the new monitoring network, by changing the location of monitoring stations, improves spatial distribution and increases forecasting accuracy of groundwater level in an area which has no data. The optimal network in the northwest of the aquifer recommends three monitoring stations. Also, a more precise survey of this area shows a wetland that has a significant impact on groundwater level in the Silakhor aquifer; therefore, recommending three monitoring stations for this area is justifiable (see Figure 7).

After designing the groundwater monitoring network using the IDW method and NSGA-II optimization algorithm, 12 optimal stations have been evaluated by the entropy method. For this purpose, a groundwater level time series in optimal stations was chosen and then it calculated the signal energy of stations in a period of time. After calculating the signals’ energy, the entropy is calculated according to Equation (1). The amount of rainfall entropy and its variations at different intervals affect the amount and diversity of the groundwater level entropy over the same time period. Comparison of the groundwater level entropy obtained from these 12 stations with rainfall entropy in the Silakhor plain indicates whether the data obtained from these stations can accurately identify these changes or these stations’ information is inadequate. To investigate this, the monthly rainfall and groundwater levels were arranged over a period of ten years. Table 2 shows the entropy value for the rainfall time series, and Table 3 shows the entropy of groundwater level value for each selected station.

Comparison of Tables 2 and 3 shows that entropy variations for each selected station are similar to the entropy variations of the rainfall time series. This indicates that the selected stations react to the entropy changes of the precipitation itself. In other words, the readings of groundwater level at these stations in different time periods are in accordance with the changes in rainfall and the results are far from real and not expected. Figure 8 shows the zoning of the study area based on the 12 optimal stations. As shown in Figure 8, the spatial distribution has been able to provide a good estimate of the groundwater level in the Silakhor plain, which results in precision and accuracy in modeling at unmatched locations as well as management plans in the context of groundwater resources management. As an example, three monitoring stations between 29 existing stations were used to estimate groundwater level by the first scenario. Table 4 shows the results of groundwater level in these points according to observation data and the first scenario. Figure 9 shows the difference between the interpolation map in the first scenario and the existing network.

The second scenario chooses optimal wells from the 29 observation wells shown in Figure 5. In the second scenario, the monitoring network was obtained from the existing observation wells in the Silakhor plain. For this purpose, 29 wells’ data are introduced into the model to determine the optimal network in accordance with the objective functions. In this scenario, the entropy of these 29 stations was calculated; then, these stations were classified as the maximum amount of entropy. After that, according to the variations of each station, the optimal network was chosen. In other words, comparison between the entropy of rainfall time series and groundwater level in 29 stations can decide on the optimal monitoring network. According to Venetsanou et al. (2016) and Abdullahi & Garba (2015) who investigated the effect of rainfall on quantitative and qualitative changes in aquifer water, the time dependency of rainfall is effective in aquifer feeding. In fact, the entropy value of the precipitation and its changes in different time intervals affect the value and variations of groundwater level entropy during the same time intervals. In this sense, the value of entropy obtained in each station in the monitoring network must be in good agreement with the time series of entropy changes. If the value of entropy obtained from the station is not proportional to the entropy variation of the time series, it indicates that the station cannot show the water balance based on the time dependence of the time series and the feeding of the aquifer.

Comparison of Tables 2 and 5 shows that entropy variations for 11 stations from the existing stations are similar to the entropy variations of the time series and, also, the entropy of these stations is higher than other stations. This indicates that the selected stations react to the entropy changes of the precipitation itself. In other words, the readings of groundwater level at these stations in different periods of time are in accordance with the changes in rainfall and the results are far from real and unexpected ones. These 11 stations are marked with stars in Table 5. The number of stations at some intervals indicates the value of entropy as zero. This value indicates that the time interval of the station in accordance with the precipitation fluctuations cannot be a precise measure of the station level in the Silakhor plain. The results of the second scenario represent a network with 11 stations among the 29 available ones. As shown in this table, 18 stations do not have the same variations according to the entropy of rainfall time series and have a few entropy variations between other stations. For this reason, these stations cannot affect the monitoring network. An investigation into the entropy of the network with 11 stations compared to the network with 29 observation stations has shown that the network with 11 stations could accommodate variations in the value of entropy according to the trend of rainfall time series. The optimal monitoring stations in the second scenario are shown in Figure 10.

According to the optimal monitoring station and data, the aquifer was zoning. Figure 11 shows the spatial interpolation map in the second scenario. As shown in Figure 11, this scenario could not predict the groundwater level as well as the first scenario on non-point in the Silakhor plain because the first scenario increased spatial monitoring accuracy by changing the location of stations. In comparison between the optimal network and existing network, it has been shown that the optimal network could reduce the cost of monitoring stations and have a similar zoning in the Silahkor plain rather than the existing monitoring network. As an example, three monitoring stations among 29 existing stations were used to estimate groundwater level by the second scenario. Table 6 shows the results of groundwater level in these points according to observation data and the second scenario. Figure 12 shows the difference between the interpolation map in the second scenario and the existing network.

The flexible nature of the first scenario as compared to the second scenario is clear in Figure 13. It is evident that optimizing the locations of the monitoring stations regardless of existing observation wells in the study area leads to the efficient collection of accurate groundwater levels.

Also, results show that taking into account the existing observation wells provides useful information for optimizing the entire monitoring network. In the second scenario, the main constraint is the fact that existing monitoring wells are used in the optimized network. This caused a rise in the RMSE in the network associated with the second scenario as compared with the RMSE associated with the first scenario (Table 7).

Considering that the implementation of water resources monitoring programs is costly and time-consuming, a method for optimizing the existing network is necessary. In the process of optimizing the groundwater monitoring network, it is important to identify the number and location of wells when faced with high volumes of data in the studied time span. A small amount of groundwater monitoring plays a special role in underground water studies and the prediction of the future behavior of groundwater flows in aquifers. Therefore, the aim of this research was to design a groundwater level monitoring network by providing an optimization method. The groundwater level monitoring network is designed by taking into account all potential points within the aquifer area. This means that the design of the monitoring network is carried out from well to well, and the search algorithm space optimizes all the potential points of the aquifer range. Proper characterization of groundwater conditions relies on well-designed groundwater monitoring networks. This study introduced the multi-objective optimization of groundwater monitoring network design with the evolutionary algorithm NSGA-II and entropy theory in two scenarios. The novel methodology was applied to the Silakhor plain aquifer, Iran. The first part of this paper created a time series of groundwater level values over the entire area of the aquifer by Kriging and stored the time series into a comprehensive database. The second part of this paper optimized the network of observation wells employing the NSGA-II and optimized groundwater levels obtained by the IDW. The objectives of the optimization problem were the minimization of the number of monitoring wells and the minimization of the RMSE.

In the first scenario, groundwater monitoring wells obtained by 12 monitoring stations have RMSE = 0.61 m. The rehabilitation of the groundwater monitoring network in the Silakhor plain has been able to improve the spatial distribution of observation wells, which was examined in the first scenario. Also, changing the location of new stations in the first scenario has increased the accuracy of the aquifer zoning at points without statistics of groundwater level. On the other hand, the new network can provide a good method for managing groundwater resources, spatial modeling of the aquifer, and subsequently, the prospecting of groundwater flow and the risk of subsidence. This scenario has no similar results to the first scenario concerning RMSE but can reduce drilling and maintenance costs.

Comparison of the alignment lines obtained from the monitoring network in different scenarios with the initial state indicates that the optimization algorithm solutions have been effective and optimal and have been able to increase the accuracy of the network based on the objective functions.

Data cannot be made publicly available; readers should contact the corresponding author for details.