Salinity assessment of groundwater for irrigation to prevent soil salinization

,

(1. Soil and Water Research Institute, Agricultural Research, Education and Extension Organization(AREEO), Karaj 3177993545, Iran; 2. Scientific Board, Soil and Water Research Institute, Agricultural Research, Education and Extension Organization(AREEO), Karaj 3177993545, Iran)

Abstract: Electrical conductivity (EC) is considered as the most important indicator for assessment of groundwater quality. Determination of suitable interpolation method for derivation of groundwater quality variables map such as EC is dependent on region conditions and existence of enough data. For determining groundwater EC, 341 groundwater samples were randomly collected from the central regions of Guilan province, paddy soils, in northern Iran. Interpolation methods including inverse distance weighting (IDW), global polynomial interpolation (GPI), local polynomial interpolation (LPI), radial basis function (RBF), ordinary kriging (OK) and empirical Bayesian Kriging (EBK) were used to generate spatial distribution of groundwater EC. The results indicate that EBK is a superior method with the least RMSE, MAE and the highest R2. The generated maps can be used to identify the regions in the studied area where groundwater could be allowed to be extracted and utilized by farmers to reduce adverse effect of the scarcity of surface water.

Key words: empirical Bayesian Kriging;geostatistics;GIS;groundwater electrical conductivity;spatial variability

Soil salinization is critically affecting our limited soil resource and deteriorating the ecosystem health. Meanwhile, it leads to soil desertification and land degradation, and results in sharp decrease in soil productivity, vegetation cover and biodiversity. Therefore land users should utilize irrigation water especially groundwater correctly. Groundwater resources are dynamic, and can continually adjust to short- and long-term in changes in climate, groundwater withdrawal and land use. The balance between charging and discharging of aquifers determines the groundwater level and its quality properties. Without a well-functioning water supply, it is difficult to imagine a productive human activity related to agriculture or livestock. Groundwater plays a fundamental role in supplying clean and safe water to competing uses for domestic, industrial and agricultural sectors such as irrigation; increasing attention is also paid to its significance for ecological integrity. However, groundwater aquifer systems are always featured by complexity, high nonlinearity, and multi-scale and random distribution as a result of the frequent interactions between surface water and groundwater as well as acute human disturbance. Thus, effective modeling techniques would be required for provi-ding efficient groundwater management strategies. Sustainable groundwater quality is important for drinking, irrigation and domestic purposes[1-3]. Groundwater resources are important water resources for agricultural uses and water drinking in Iran and many other countries which have climates similar to Iran′s[4]. Due to less pollution of groundwater sources, even in areas where there is sufficient surface water resources, groundwater is still used to irrigate crops. Toxic agents can be produced from fertilizers, pesticides and se-wage[4-5]. Therefore, during recent years, increasing pollution and losing of water resources have changed exploitation policy of water and soil resources[6].

Prior to the design of groundwater quality networks, it is essential to investigate the spatial structure of the groundwater quality variables monitored such as electrical conductivity. Generally, the objective of monitoring introduces these variables. For example, greater effort is required for monitoring groundwater used for domestic (municipal) purposes than for agricultural use. The aim of characterizing the spatial structure of the variables is to assess their monitoring, which is needed to assess the cost of monitoring, and to give a clear picture about their spatial variability or structure. The spatial structure of the groundwater quality variables can produce, for example, contour maps of the variable means. These maps can be used for predicting and signifying pollution areas. Accordingly, protection measures, and management and planning decisions can be made to minimize the deterioration in the polluted areas[7-9].

Different statistical and geostatistical approaches have been used in the past to estimate the spatial distribution of groundwater electrical conductivity or sali-nity[10-11]. Classical statistics could not make out the spatial allocation of soil and water properties at the unsampled locations[12]. Geostatistics is an efficient method for studying the spatial allocation of soil cha-racteristics and their inconsistency. It can reduce the variance of assessment error and execution costs[13-14]. ZEHTABIAN, et al[7]used two techniques including Kriging and weighted moving average (WMA) techniques for presenting spatial variation of groundwater properties such as electrical conductivity. Finally, comparison of the results using statistical techniques shows that Kriging technique performs better than WMA technique. BARAM, et al[15]demonstrated that analyzing vadose zone and groundwater data by spatial statistical analysis methods could significantly contribute to understanding the relationships between groundwater contaminating sources, and to assessing appropriate remediation steps. BODRUD-DOZA, et al[1]used ordinary Kriging interpolation method for taking initial decision of spatial distribution of groundwater quality parame-ters. Their results represented that the ordinary Kriging technique was capable of predicting spatial variability more accurately for the studied area with suitable semivariagram model. They reported that outcomes of the study would provide insights for decision-makers to take proper measures for groundwater quality management in central Bangladesh. HUSSAIN, et al[16]reported that empirical Bayesian Kriging (EBK) was most suitable for spatial prediction of total dissolved solids (TSD) in drinking water. MIRZAEI, et al[17]reported that EBK model was the best of all the geostatistical models including ordinary Kriging (OK) and inverse distance weighting (IDW) for estimation of groundwater contamination. GUNARATHNA, et al[2]used IDW, radial basis function (RBF), Kriging and EBK for interpolation of groundwater salinity. They reported that EBK method had the least root mean square error (RMSE) value for the spatial variation of groundwater electrical conductivity. ABU-ALNAEEM, et al[9]used an integrated statistical, geostatistical and hydrogeochemical approach for assessing groundwater salinity and quality in Gaza coastal aquifer, Gaza Strip, Palestine. JASROTIA, et al[18]prepared the spatial distribution map of various physical-chemical parameters using the geographic information system (GIS) and delinea-tion of groundwater quality zones in the Lesser Himalayan region.

Groundwater monitoring of can provide fundamen-tal information for sustainable water resource management for irrigation to prevent soil salinization. The goals of groundwater monitoring can be ambient resource condition, compliance, risk detection, research monitoring and especially water salinity, or a combination of these goals. Land and water management practices should be developed according to results of continuous monitoring of water table depth and groundwater qualities. In irrigated areas, monitoring of wells are commonly used to evaluate spatial and temporal changes in water table level and groundwater quality such as salinity. Groundwater quality variables mapping is an important tool for groundwater management and risk assessment. In recent years, agricultural fields in Guilan province in northern Iran was dominantly irrigated by groundwater, and in recent decades, exploitation of water and soil resources has been generally changed by excavation of many deep and mid-deep wells. Therefore, sustainable management of water and soil resources requires to be informed of groundwater quality changes. The present study was carried out with objective to evaluate the accuracy of different interpolation methods including inverse distance weighting (IDW), global polynomial interpolation (GPI), local polynomial interpolation (LPI), radial basis function (RBF), ordinary Kriging (OK) and empirical Baye-sian Kriging (EBK) for predicting groundwater electri-cal conductivity parameters in central areas of Guilan province in northern Iran for suitable management of irrigation water to prevent soil salinization. The results of our study could be used by the stakeholders to perform a more complete and explanatory plan for the optimum management of soil and water resources.

1 Materials and methods

1.1 Studied area

The studied area is located between 49°31′-49°45′ E longitude and 37°7′-37°27′ N latitude in Guilan province of northern Iran bordering to the Caspian Sea. The climate of the region is humid with the mean annual precipitation of 1 293.6 mm. The mean annual temperature of the region is 15.8 ℃. The mean humidity is 75% and the annual evapotranspiration is 850 mm. The soil moisture and temperature regimes of the region by means of Iran regimes maps are Udic or Aquic and Thermic, respectively. The major geological formations are composed of thick sedimentary and metamorphic rocks of Tertiary and Quaternary periods. The coastal plain lying between the Alborz mountain ranges and the Caspian Sea is composed of marine, river and aeolian deposits of varying thicknesses. The physiographical units of the region from south to north direction are upper plateaus, river alluvial plains, river bank, lowlands and coastal lands, respectively. This studied area is 40 000 hm2. A total of 341 samples were collected from groundwater randomly for determining electrical conductivity. Electrical conductivity was mea-sured using a meter and probe as well. The probe consists of two metal electrodes spaced one centimeter apart (thus the unit of measurement is microSeimens or milliSeimens per centimeter). A constant voltage was applied across the electrodes resulting in an electrical current flowing through the aqueous sample. Fig.1 shows the distribution of sampling points in the studied area.

Fig.1 Location of studied area and sampling points

1.2 Interpolation methods

In the present studies, deterministic (i.e., create surfaces from measured points) and geostatistical (i.e., utilize the statistical properties of the measured points) interpolation techniques were used. In this study, a variety of deterministic interpolation techniques, including those based on the extent of simila-rity IDW, LPI, GPI, degree of smoothing (RBF) and geostatistical interpolation, namely OK and EBK were used to generate the spatial distribution of groundwater electrical conductivity.

1.2.1 Inverse distance weighting

The inverse distance weighting(IDW) is one of the mostly applied deterministic interpolation techni-ques in the field of soil science. IDW estimations were made based on the nearby known locations. The weights assigned to the interpolating point are the inverse of its distance from the interpolation point. Consequently, close points are made up to have greater weights (so, more impact) than distant points and vice versa. The known sample points are implicit to be self-governing from each other[14].

1.2.2 Radial basis function

Radial basis function (RBF) methods are a series of exact interpolation techniques, that is, the surface must go through each measured sample value. RBF me-thods are a form of artificial neural networks. RBFs are like a rubber membrane which is fitted to each measured data point while minimizing the total curvature of the surface must pass through each sampled point[14].

1.2.3 Global polynomial interpolation

Global polynomial fits a polynomial formula to the sample points. Conceptually, global polynomial posi-tions a plane between the sample points. The unknown height is then determined from the value on the plane which corresponds to the prediction location. The goal for global polynomial is to minimize errors. Global polynomial interpolation (GPI) fits a smooth surface defined by a mathematical function (a polynomial) to the input sample points. The global polynomial surface changes gradually and captures coarse-scale pattern from the data. Conceptually, GPI is like taking a piece of paper and fitting it between the raised points (raised to the height of value).

1.2.4 Local polynomial interpolation

Local polynomial fits many smaller overlapping planes to the sample points, and then uses the center of each plane as the prediction for each location in the studied area. Local polynomial interpolation(LPI) creates a surface from many different polynomial formulas, and each is optimized for a specified neighbor-hood; the neighborhood shape, maximum and minimum number of points, and a sector configuration can be specified; the sample points in a neighborhood can be weighted by their distance from the prediction location.

1.2.5 Ordinary Kriging

Kriging is one of the most popular and robust interpolation techniques. It integrates both the spatial correlation and the dependence in the prediction of a known variable. Estimations of nearly all spatial interpolation methods can be represented as weighted averages of the sampled data[1]. The presence of a spatial structure where observations close to each other are more alike than those which are far apart (spatial autocorrelation) is a prerequisite for the application of geostatistics. The experimental variogram measures the average degree of dissimilarity between non-sampled values and nearby data value, and thus can depict autocorrelation at various distances[13].

Variogram plots (experimental variograms) were acquired by calculating variogram at different lags. Gaussian model was selected in order to model experimental variogram and to acquire information about the spatial structure as well as the input parameters for Kriging estimation.

From the analysis of the experimental variogram, a suitable model (e.g. Gaussian, spherical, and exponential) is then fitted usually by weighted least squares and the parameters (e.g. range, nugget and sill) are used in the Kriging procedure. The ratio of nugget effect to sill can be considered for evaluation of spatial dependence of data. When this ratio is smaller than 0.25, the concerned parameter has a strong spatial dependence. The spatial dependence between 0.25 and 0.75 is middle; when it is greater than 0.75, the spatial dependence is weak[14].

1.2.6 Empirical Bayesian Kriging

Empirical Bayesian Kriging (EBK) automates the most difficult aspects through a process of subsetting and simulations. EBK process implicitly assumes that the estimated semivariogram is the true semivariogram for the interpolation region and a linear prediction that incorporates variable spatial damping. The result can be obtained through a robust non-stationary algorithm for spatial interpolating geophysical corrections. This algorithm extends local trends if data coverage is good while it allows for bending to a priori background mean if data coverage is poor[14].

1.3 Data preprocessing for use in prediction methods

The first step for using geostatistic methods is to study the existence of spatial correlation between data by variogram analysis. The condition of this analysis is that data must be normal. One of the evaluation methods for nominating normality of data is the usage of skewness coefficient. When skewness coefficient is lower than 0.5, there is no need to convert data, however, if this coefficient is between 0.5 and 1, and more than 1 for normalizing data, square root and lo-garithm must be used, respectively. In order to know whether the data were normal or not, Kolmogorov-Smirnov test was used.

1.4 Validation analysis of methods

Three different types of standard statistical perfor-mance evaluation criteria were used to control the accuracy of the prediction capacity of the models developed. These are root mean square error (RMSE), the coefficient of determination (R2) and mean absolute error (MAE). Performance evaluation criteria used in the current study can be calculated using the following equations:

(1)

(2)

(3)

All statistical calculations were performed using Microsoft Excel 2007 and SPSS 24.0. Geostatistical analysis and generation of prediction maps of water electrical conductivity were carried out with ArcGIS 10.3.1 software (ESRI, Redlands, CA, USA).

2 Results and discussion

2.1 Descriptive statistics

Some statistical characteristics such as mean, standard deviation, minimum, maximum, coefficient of variation, skewness and kurtosis are presented in Tab.1 for the variables of electrical conductivityσand its square root transformσ*. Since skewness of data is more than 0.5, data do not obey normal distribution. Therefore, square root transform was used for data normalization and its results are mentioned in Tab.1. These results show that data are normalized because skewness is smaller than 0.5. The Kolmogorov-Smirnov test also demonstrates that the transformed data obey normal distribution. For electrical conductivity parameterσ, an analysis trend was made, and it was determined that there was no global trend.

Tab.1 Statistics of data used in research

2.2 Interpretation of prediction methods

Empirical Bayesian Kriging was performed for mapping groundwater salinity. EBK variogram is presented in Fig.2. The fitted model is K-Bessel detrended. It had suitable results in the prediction of groundwater salinity compared to other models. In this method, transfor-mation type was log empirical, with a smoothing factor of 0.2 and a radius value of 2 486.16 m.

In the semivariograms of the examined water parameters, the blue crosses represent the semivariance obtained from empirical transformation, the bold red lines represent the intermediate distribution, the faded red plots (dotted) display the 25th and 75th percentile respectively, and the densely-packed blue lines indicate a series of semivariogram passing through the zone

Ordinary Kriging technique was performed for interpolation of groundwater electrical conductivity data. As seen in the variogram results (Fig.3 and Tab.2), the most appropriate model fitted to groundwater electrical conductivity is Gaussian. The model fitted to empirical semivariogram of groundwater electrical conductivity has characteristics including active lag distance 10 000, lag class distance interval (uniform interval) 800, offset tolerance degree 22.5°, and neighbors to include 20 including at least 15 and sector type perpendicular tow line, i.e. add symbol.

Fig.3 Experimental semivariogram of groundwater salinity and its fitted model by ordinary Kriging

Since the model fitted to empirical semivariogram has highR2and lowRSS, the model fitted to empirical semivariogram is the most appropriate among other models, i.e., IDW, RBF, LPI and GPI models. The ratio of nugget effect to sill is smaller than 0.25 for va-riable; therefore the fitted Gaussian model has a strong spatial dependence[14]. The prediction maps of ground-water salinity by OK and EBK are presented in Fig.4.

IDW method was performed on electrical conductivity data with power values from 1 to 5. This model had neighbors to include 14 and circular sector type. In this method, the best result was derived from power value equaling to 1 which was more accurate than other power values. The fitted RBF model for groundwater electrical conductivity data has characteristics including kernel function spline with tension, parameter 0.09, and neighbors to include 16 and circular sector type. The prediction maps by these methods are presented in Fig.5.

Tab.2 Properties of used model in ordinary Kriging method

Fig.4 Spatial distribution map of groundwater electrical conductivity using OK and EBK methods

The characteristics of used LPI technique in this research include optimized weight distance with a weight of 16 724.409, and neighbors to include 341, including at least 10 and sector type multiple. In the prediction of electrical conductivity data, results obtained by this method with power value equaling to 4 are more accurate than those by IDW. GPI prediction model for groundwater electrical conductivity data had the best export with power value equaling to 6 com-pared with others. The spatial distribution maps of groundwater electrical conductivity using these methods are presented in Fig.6.

2.3 Comparison among different interpolation techniques

IDW, RBF, LPI, GPI, OK and EBK were used to estimate groundwater electrical conductivity. After evaluating different models, it was concluded that the Gaussian model had the most accurate prediction and therefore, it was selected as the best model fitted to the data in OK method. In EBK method, the best fitted model was K-Bessel detrended. The summary statistics by geostatistic method showed that empirical Bayesian Kriging with K-Bessel detrended model provided much better estimation results for electrical conductivity than other methods (Tab.3), and this result was in agreement with the findings of HUSSAIN, et al[16], MIRZAEI, et al[17]and GUNARATHNA, et al[2]. Kriging is a widely used method of geostatistical interpolation which assumes that no regional trend exists in the data. By comparison among different methods, RBF method had better result than IDW to simulate groundwater electrical conductivity variable. Results obtained by IDW with power value equaling to 1 are more accurate than those by LPI and GPI.

Fig.5 Spatial distribution map of groundwater electrical conductivity using IDW and RBF methods

Fig.6 Spatial distribution map of groundwater electrical conductivity using LPI and GPI methods

Results of the current study showed a strong spatial dependence of the variable data, but the most appropriate results based on statistical comparisons showed high capability of Kriging techniques especially EBK because its statistical characteristics such as coefficient of determination, root mean square error and mean absolutely error were better than those of other prediction methods (Tab.3). Generally, our results were similar to the findings of ZEHTABIAN, et al[7], BHUNIA, et al[14]and BARAM, et al[15].

The validation and sufficiency of the developed model variogram can be tested via a technique called cross validation. In this method, estimation is obtained by leaving one sample out and using the remaining data. This test allows evaluating the goodness of fitting of the variogram model, the appropriateness of neighborhood and type of Kriging used. The interpolation values are compared with the real values and then the least square error models are selected for regional estimation[13]. Cross validation results of predicted electri-cal conductivityσpand measured electrical conducti-vityσmby EBK method are presented in Fig.7. This scatter plot shows well evaluation of estimation by EBK method.

Spatial distribution maps of groundwater electrical conductivity estimation values by different interpolation methods are shown in Fig.4-6. Groundwater salinity is higher in the northern of the studied area, near the Caspian Sea in coastal land and lowlands. High electrical conductivity of coastal land groundwater is due to the seepage effect of sea water while in lowlands it is because of surface drains accumulation from adjacent areas and irregular use of chemical fertilizers.

Tab.3 Evaluation results of prediction by different methods

Fig.7 Scatter plot of measured versus predicted electrical conductivity using EBK method

Therefore,it is suggested that effective measures be taken to prevent the groundwater electrical conductivity from increasing in lowlands and coastal land of the studied area. Drainage is necessary for water removal in case of rising water table and for removing the salts in case of increasing the salinization in root zone. Land drainage is one of the key inputs to get better yields per unit of farmland. Therefore, salinization and drainage problems of irrigated agriculture could be considerably affected via suitable management in the studied area. Thus, if effective procedures were not carried out for the control of electrical conductivity increasing, these lands will be lost in future for farmers cannot cultivate in those lands because of their low yields and high costs.

2.4 Analysis of groundwater salinity and optimal management for preventing soil salinization

Improper and excessive use of irrigation water is one of the major factors aggravating the groundwater salinity. CHAUDHURI, et al[19]reported that irrigation return flow was a major mechanism of solute enrich-ment of groundwater systems in agricultural regions in the Ogallala aquifer in the United States, and the eva-porative enrichment of salts in the upper part of soil profile and subsequent leaching of salts with irrigation water was the main reason for salt enrichment of shallow groundwater systems in Ogallala aquifer. Moreover, chemical species causing high electrical conductivity and their potential sources should be identified to mitigate their further damage on groundwater and soils. In the above-mentioned area of the United States, irrigation water return flow has a diluting effect on electrical conductivity of groundwater so that electrical conductivity is generally reduced following irrigation seasons in lowlands; these results are in line with the findings of previous studies. The United State Salinity Laboratory classified groundwater into four grades by electrical conductivity values: groundwater is excellent if electrical conductivity valueσ<0.25 dS/m, good ifσranges from 0.250 to 0.750 dS/m, fair ifσis between 0.75 and 2.25 dS/m, and poor ifσ>2.25 dS/m. Our results exhibit that the groundwater electrical conductivity of the majority of the studied area falls into good category, except some locations in the north of the studied area. Interactions of irrigation water with natural processes should be recognized in the analysis of groundwater salt enrichments. Under intensive irrigated agricultural production, a considerable amount of salts may move beyond the root zone, degrading the groundwater. Plants uptake nearly pure water and leave salts behind; the salts are then transported to ground-water by percolating water. During transportation, the salt-rich percolating water interacts with soil and rock constituents, and thus releases chemical species, further rising the salt concentration of groundwater. The above-mentioned factors must be carefully controlled to minimize the problems.

Identification of spatial and temporal pattern of groundwater salt concentration and groundwater salinity is an important step in setting appropriate alternative management practices to protect land and soils against degradation. The data obtained in this study may help mitigate soil and groundwater degradation caused by climate change and human influence. Particular attention should be paid on the locations with high groundwater electrical conductivity. The results may have important implications for other countries with similar climate, topography and soil conditions. It is important to define the effect of irrigation and agricultural water return flow in combination with chemical fertilizers on the quality and quantity of groundwater in the studied area. In addition, the impact of sea water intrusion from the Caspian Sea via the bed of Anzali Lagoon should be frequently monitored to avoid soil and groundwater degradation due to increasing salinity in coastal plain and lowlands in the studied area. ABUALNAEEM, et al[9]argued that sea water intrusion in coastal parts, vertical and lateral mixing of water, and anthropogenic inputs are responsible for groundwater salinization. For the control of water quality, rivers should be monitored in regard with climate change to avoid irreversible impact of irrigation on soils and groundwater as already observed in many locations of Iran. Degraded quality of groundwater due to increasing salinization and salt ions of fertilizers was apparent from our results. A combination of natural and anthropogenic processes has caused salinization in shallow ground-water in the studied area. In some instances, natural processes were triggered by anthropogenic sources such as fertilizers, irrigation and domestic waste disposals. Our results showed that an increase in the concentration of fertilizer salt ions and high salinity of groundwater were growing concern in the studied region. In addition to aquifer quality, shallow groundwater contamination in the studied region should be considered in developing and implementing strategies for rural development. These conclusions should also be carefully considered in order to prevent soil salinization.

3 Conclusions

In this study, geostatistical methods are used to analyze spatial dependence and changes of groundwater electrical conductivity. The following conclusions are obtained:

1) Research of spatial variation of groundwater electrical conductivity is necessary for optimal management of groundwater resources.

2) According to evaluation criteria, the accuracy of geostatistical methods in the estimation of ground-water electrical conductivity is very high.

3) Empirical Bayesian Kriging can be used as a suitable tool to estimate the salinity of groundwater in the studied area with data restriction.

4) Empirical Bayesian Kriging method has higher accuracy in estimation of groundwater electrical conductivity than other methods. The groundwater electri-cal conductivity is high in coastal land and lowlands compared to other areas, therefore it should be noticed and effective measures are suggested to be taken to prevent electrical conductivity value from increasing.

5) Generally, geostatistics models are suitable for estimation of groundwater quality variable and geostatistical interpolation methods such as ordinary Kriging and empirical Bayesian Kriging perform better than deterministic interpolation methods, for example, RBF, IDW, LPI and GPI in mapping the spatial distribution of groundwater salinity.

4 Acknowledgement

The authors thank the editor and the anonymous reviewers for their suggestions for improving the quality of the manuscript. The authors wish to express their sincere thanks to all members of the Soil Science Laboratory of Soil and Water Research Institute (SWRI). We are grateful for the support from Water and Soil Deputy of Jihad-e-Agriculture Ministry.