Introduction

Geographical studies are basic tenets for environmental planning. Understanding an area’s attributes can make it easy for reasonable measures. Iran is located in an arid and semiarid area of the world which is immediately affected by random rainfalls. Based on Climate Atlas during 1961–1990, more than two-thirds of the country had less than 300 mm rainfall annually.¹ Currently, one of the main techniques to achieve a sustainable development in planning and also preventing from damages such as flood, erosion, and frost is through measuring meteorological parameters. This aim by traditional method requires a network of rainfall station, but topographical situations are costly which may deter countries to make such system in a large number. Hamadan province as the subject of this study has wide watersheds with different hydrological conditions. Furthermore, due to the few number of hydrometric and rainmetric stations in the region, statistics and information on rainfall and runoff cannot precisely show place changing all around the province. That is why accurate information on interpolation is highly important. Interpolation is an appropriate approach that has a widespread use of generalizing spot data to region data.² The main difference between this method and the classical statistics method is that in classical statistics, samples drawn from a population are independent, which means each sample does not give information about the next sample, whereas geostatistical methods show varied accuracy of estimating a neighboring point. Thereby, using these highlights their vitality in geographical information systems and selecting appropriate techniques. In 1990, geostatistical methods were used to optimize monitoring networks in Morocco. The results indicated the superiority of geostatistical methods in comparison with classical statistics.³ In California, a research was done using 2 hydrological parameters and a statistical method (co-kriging). The result of this research markedly specified more accurate estimation than traditional techniques.⁴ McKenna used geostatistical techniques that were based on variogram to enhance the accuracy of the estimation and simulation of data on groundwater flow in Colorado. The results of this examination showed that the use of geostatistical methods increases the accuracy of data; similarly, these techniques require a few number of samples.⁵ Misaghi simultaneously used the geostatistics and artificial neural network algorithm to extract the spatial distribution of rainfall for the area in Khuzestan province. The results illustrated that among all the algorithms, geostatistical approach has a better result. Also, error estimates for neural network models were greater than geostatistical methods.⁶ In a separate study, Carrera-Hernandez used various forms of geostatistical method to analyze daily climatic data from approximately 200 stations located in the Basin of Mexico for the months of June in 1978 and 1985. The results of cross-validation showed that the interpolation of daily events was improved by the use of elevation as a secondary variable even when that variable showed a low correlation.⁷ Rahimi used geostatistical methods to estimate rainfall in an arid and semiarid area in the southeast of Iran. The results demonstrated that for annual rainfall methods such as thin plate smoothing splines (TPSS) with altitude auxiliary variables and kriging were qualified approaches, whereas for monthly rainfall TPSS method with a power of 2 was the best method.⁸ In another study, different methods to reconstruct the statistical gaps of rainfall in the Central Alborz region were applied. The methods in the study included linear regression, normal ratio, inverse distance weighting, and geostatistical methods. Meanwhile, among the total number of stations, 18 stations with full statistics and no gap were selected. The results indicated that the geostatistical methods are the most appropriate ways for correlation coefficient.⁹ Hooshmand used statistical methods such as kriging and co-kriging for simulating water quality variables and concluded that co-kriging is a fully qualified geostatistical method for the simulation of groundwater quality variables.¹⁰ In another study, geostatistical techniques were used to assess the effect of incorporating topographical data with geostatistical interpolation for monthly rainfall and temperature in Ping Basin, Thailand. It was found that co-kriging models of mean monthly rainfall and mean monthly temperature are more effective than other interpolations.¹¹

Methods

The study area is in Hamadan province which is one of the mountainous areas located in the west between latitudes 33°33′ to 35′48′N and longitudes 47°45′ to 49′36′E. It has an area of about 19 545/82 km². Alvand mountain with an altitude of 3574 m above sea level caused mountainous climate with cold winters and mild summers. Its average rainfall and temperature are 317 mm and 30°C, respectively. To determine the amount of rainfall, data from 35 synoptic stations and rain measurements are collected, which have no defect and are shared over 30 years from 1982 to 2012. The location of Hamadan province in Iran and characteristics of each station are illustrated in Figure 1 and Table 1, respectively.

Figure 1.

Location of the study area.

Table 1.

Stations’ characteristics.

Spatial prediction technique

Kriging

Kriging methods rely on the notion of autocorrelation. There is difference between classical statistics and geostatistics technique. Over classical statistics, observations are assumed to be independent which means there is no correlation among observation samples. In geostatistics, autocorrelation is known as a distance function. The information on spatial locations allows researchers compute distances among observation samples and model autocorrelation as a function of distance. This is related to a basic principle of geography issues: Observation samples closer together tend to be more similar than those that are farther apart.^12,13 The experimental variogram measures the average degree of dissimilarity between unsampled values and a nearby data value and thus can depict autocorrelation.¹⁴ The value of the experimental variogram for a separation distance of h (referred to as the lag) is half the average squared difference between the value at z(x_i) and the value at z(x_i + h)¹⁵:

where N(h) is the data pair number within a given class of distance and direction. If the values at z(x_i) and z(x_i + h) are auto-correlated, the result of equation (1) will be small, relative to an uncorrelated pair of points. From an analysis of the experimental variogram, a suitable model is fitted, usually by weighted least squares, and the parameters (eg, range, nugget, and sill) are then used in the kriging procedure. In this study, ordinary kriging and simple kriging are used with their models.

Co-kriging

Co-kriging is defined as a multivariate version of kriging that minimizes the variance of estimation error by considering the spatial correlation between the variables of interest and secondary variables. The cross-semivariogram is used to quantify cross-spatial autocovariance between the original variable and the covariate.¹⁶ The cross-semivariance is computed by equation (2):

where ${\gamma }_{uv}h$ is the cross-semivariance between u and v variables, $Z_u(x)$ is the primary variable, and $Z_v(x)$ is the secondary variable.

Radial basis functions

Radial basis function (RBF) methods are “a series of exact interpolation techniques; that is, the surface must go through each measured sample value. Each basis function has a different shape and results in a slightly different interpolation surface.” Radial basis function showed that radial basis function neural network gives the same results as a standard geostatistical technique: prediction mapping of the regression function and mapping of the prediction variance. Radial basis function neural network is expressed by equation (3):

where lgr; is a zero-mean noise and h is the RBF.

Radial basis functions are a class of functions, the key feature of which is that the distance from the center determines their response. Most commonly used RBF is the Gaussian RBF. Gaussian model provides an example of local RBF that gives significant response only in a neighborhood of its center.^171819–20

Inverse distance weighting

Inverse distance weighting (IDW) is “a deterministic, nonlinear interpolation technique that uses a weighted average of the attribute (ie, phenomenon) values from nearby sample points to estimate the magnitude of that attribute at nonsampled locations.” The weight a particular point is assigned in the averaging calculation depends on the sampled point’s distance to the nonsampled location. These weights are controlled on the bases of power of 10. Often with an increase in power, the effect of the points that are farther diminishes. Lesser power distributes the weights more uniformly between neighboring points. We should keep in mind that in this method the distance between the points counts, so the points of equal distance have equal weights.^21,22 In this method, the weight factor is calculated by equation (4).

where ${\mathrm{lgr;}}_i$ is the weight of point, $D_i$ is the distance between point, $i$ is the unknown point, and α is the power 10 of weight.

Global polynomial interpolation

Global polynomial interpolation (GPI) is “a quick deterministic interpolator that is smooth. There are very few decisions to make regarding model parameters. It is best used for surfaces that change slowly and gradually.” However, there is no assessment of prediction errors and it may be too smooth. Locations at the edge of the data can have a large effect on the surface. There are no assumptions required of the data.

Local polynomial interpolation

Local polynomial interpolation (LPI) is “a moderately quick interpolator that is smooth. It is more flexible than the global polynomial method, but there are more parameter decisions.” The method provides prediction, prediction standard error, and condition number surfaces that are comparable with ordinary kriging with measurement errors.

Comparison of different methods

To evaluate interpolation methods, statistical indices of mean absolute relative error (MARE) and correlation coefficient are used. The lowest MARE indicates the greatest reliable prediction. The MARE was computed according to equation (5):

where Z*(X_i) is the predicted value at point X_i, Z(X_i) is the observed value at point X_i, and n is the number of samples,

Results and Discussion

In this study, various interpolation methods to estimate rainfall distribution were used, the results of which are as follows. At the first step, we examined the normality of data by Kolmogorov-Smirnov test in SPSS. The result showed an acceptable normality of the data. Then, among the used methods, simple co-kriging illustrated better result for rainfall distribution in the mentioned area. According to the correlation coefficient and MARE, the best models are chosen from geostatistical methods, namely, simple co-kriging (exponential), ordinary co-kriging (exponential), local polynomial interpolation (model with power 4), global position system (model with power 1), RBF (inverse multiquadric), ordinary kriging (Gaussian), simple kriging (hole effect), and inverse distance weighting (model with power 1). Tables 2 and 3 show results of interpolation methods. The results of this study are confirmed by researchers including Safari,²³ Golmohammadi,²⁴ and Hooshmand,¹⁰ who reached the same consequence in their studies. These specialists declared that co-kriging is a qualified technique for their study areas; it should be mentioned that these research works were done in different subject areas of Iran under similar weather conditions. In contrast, Bajjali²⁵ reached a different result of his study when he used geostatistical methods to examine artificial recharge dams and its effect on the quality of groundwater and found kriging and IDW techniques the best when evaluating the performance of the artificial dams in coastal areas, but in this study, kriging and IDW techniques showed low accuracy; this can be related to its moderate weather that has privileged 317 mm rainfall annually in comparison with coastal areas that seasonal rains have played an effective role in these areas. Figures 2 and 3 show the spatial distribution of precipitation using the best models of the mentioned techniques. Turning to details based on co-kriging’s figure, it can be realized that the north and the northeast of the region are the area with low rainfall during the study period which privileged an average 1812 m height above sea level which include parts such as Famenin, Nojeh, Razan, Kitu, Omarabad, Zahtaran, Gharakand, Aab Baric, Subashie, and Shirin SU. However, the west, south, southeast, and southwest of the area showed much more rain, which contain geographical areas such as Taric Dareh Bala, Cheshmeh-Ali, Khoram Rood, Farsban, Valashgerd, Vandarabad, Asadabad, Hamadan, Nahavand, and Tuyserkan with an average 1866 m height above sea level. It can be obviously realized that there is a correct correlation between the accuracy of co-kriging technique and elevation changes. Carrera-Hernandez and Gaskin⁷ declared that the interpolation of daily events was improved by the use of elevation as a secondary variable even when that variable showed a low correlation.

Table 2.

Results of various methods of interpolation.

Table 3.

Results of various methods of interpolation.

Figure 2.

Spatial distribution of annual precipitation using kriging methods: (A) ordinary kriging, (B) simple kriging, (C) simple co-kriging and (D)ordinary co-kriging.

Figure 3.

Spatial interpolation of annual precipitation using (A) inverse distance weighting, (B) global polynomial interpolation, (C) local polynomial interpolation and (D) radial basis function.

Conclusions

To achieve an accurate and logical result in every research, it is necessary to examine all methods to reach a suitable one. In such situations, planning gives the best result to avoid possible damage and make optimal use. The main goal of this study is assessment of rain distribution and also the effect of considering elevation as a secondary variable to interpolate rainfall. Hamadan province with its wide area and random rainfall, especially in fall and winter, and landslide-prone slopes and lack of knowledge about the amount of rainfall all over the region makes it a vulnerable zone. To solve regional problems and make correct decisions, the interpolation of rainfall using GIS in the area was applied. Thus, among the aforementioned methods, simple co-kriging (exponential) is recommended with the highest correlation coefficient of .75 and the lowest amount of MARE of .12. Furthermore, according to the simple co-kriging’s figure, it can be seen that areas such as Taric Dareh Bala and Cheshmeh-Ali with an average of over 2000 m height above sea level which are located in the southeast enjoyed 542.2 and 509.9 mm of rainfall, respectively. Besides, Rasoolabad with the highest average rainfall of 658.74 mm and 2031 m height above sea level significantly proved the impact of altitude changes on the amount of rainfall. Similarly, areas such as Qahavand and Famenin which are located in east and northeast of the study area by 1715 and 1727 m height above sea level, respectively, enjoyed rainfall of about 240 mm annually. As a consequence, it can be concluded that the rain amount is immediately influenced by elevation, and also the result of simple co-kriging (exponential) produced more realistic smooth interpolated surface maps. That is to say, the result of co-kriging has a correct correlation with elevation changes unlike other geostatistical techniques such as LPI, GPI, RBF, kriging and IDW which showed low accuracy in their assessments. Therefore, they are not considered accurate interpolators.