Measurement and modeling of particulate matter concentrations: Applying spatial analysis and regression techniques to assess air quality

Graphical abstract

combined techniques can overcome the limitations of all other aerosol monitoring products. The Haz-Dust EPAM-5000 uses the principle of near-forward light scattering of an infrared radiation to immediately and continuously measure the concentration of airborne dust particles in mg/m 3 . Sampling stations were selected according to the study manual for the European Study of Cohorts for Air Pollution Effects [3]. The monitoring sites were located at 25 m of a traffic light or a traffic junction, with a sampling height at least 1.5 m above the ground. Furthermore, monitoring sites were not located within 25 m of locations where smokers are allowed to smoke and/or gather (for example close to the entrance of restaurants, hospitals, schools or other public buildings), because smokers usually congregate close to the entrance of these buildings. Instrument calibration was done according to the procedure adopted by Environmental Devices Corporation. A GPS (Global Positioning System) instrument with UTM system (Model eTrex Vista) was used for geographical positioning (X: longitude and Y: latitude). The names of the sampling sites and their related geographical coordinates in Sabzevar are listed in Table 1. Concentrations reported in Table 1 were measured in a 9-month study period.

Interpolation methods
Interpolation is a procedure to predict the value of attributes at non-sampled sites from measurements made at point locations within the same area. Interpolation is used to convert data from observation points to continuous fields to compare spatial patterns sampled by these measurements with spatial patterns of other spatial entities. The rationale behind spatial interpolation is the very common observation that, on average, values at points close together in space are more likely to be similar than points further apart. Among spatial interpolation methods, one can find RBF, IDW, and Kriging techniques [4]. RBF methods predict values that can vary above the maximum or below the minimum of the measured values. For all RBF methods, a parameter controls the smoothness of the resulting surface. By using radial basis functions, dealing with higher dimensional problems in a similar way as is possible dealing with two-and three-dimensional problems. Splines (RBF) are interpolators fitting a function for sampled points. The algorithm uses a linear combination of n functions, one for each known point as demonstrated in the following equation: where '(r) represent the interpolation function, S i À S 0 the Euclidean distance r between an unknown point S 0 and an observed one S i , while v i , with i = 1,2, . . . ,n + 1, are weight. Weights are assigned according to the distance of known points, under the constraint that, in their locations, the function must give the measured value. Radial Basis Functions (RBF) is a family of five deterministic exact interpolation techniques [5][6][7] as below: Thin-plate Spline function: '(r) = (s.r) 2 ln (s.r) Multi-quadric function: 'ðrÞ r 2 þ s 2 Inverse Multi-quadric function: Completely regularized Spline function: 'ðrÞ ¼ À X 1 n¼1 À1 ð Þ n :r 2n n!n ¼ ln Spline with tension function: where r = distance between the point and the sample, s = tension factor, E 1 = exponential integral function, C E = constant of Eulero (0.577215), and K 0 = modified Bessel function. Splines functions are slightly different, each one has a different smoothing parameter depending on the s parameter. In every method, the higher the value of s, the higher the gradualness of the variation, except for the inverse multi-quadric where the opposite condition is true. The regularized Spline creates a smooth, gradually changing surface. The regularizing parameter is in fact employed to achieve a smoother solution: e.g. a small value results in a close approximation of the data, while a large one results in a smoother solution [8]. IDW interpolates all values of the points within the sample range as averaging tool and gives better interpolation estimates when the minimum and maximum values of the surface are represented by sample data points. The concentration of point will have heavier weight if it is proximal to the required point and vice versa. Here, weight is an inverse function of the distance, as demonstrated in the following equation: where Z j is the concentration at the j th point, W i is the weight of observed i th point, d ji is the distance from the i th point to the j th point, p is the power and n is total number of points [9]. The IDW functions used in this study are power 1 and power 2.
Recognizing that the spatial variation of any continuous attribute is often too irregular to be modelled by a simple, smooth mathematical function, Kriging is a wide family of interpolation methods using geostatistics. Geostatistical methods for interpolation rely on the assumption of spatial autocorrelation. This suggests that the distance and direction between sample points are the major factors governing the estimated values at unknown points. With spatial autocorrelation, the function used for data fitting provides better estimates. Kriging methods include Ordinary, Simple, Universal, Indicator, Probability, Disjunctive Kriging, and Co-Kriging, which all rely on the concept of autocorrelation. The Ordinary Kriging Eq. (8) represents a weighted sum of the data, which is: where Z S 0 ð Þ is the predicted value at location S 0 , N is the number of samples used for the estimation, Z S i ð Þ is the measured value at location i, l i is the unknown weight for the measured value at i th location determined using the fitted variogram [4]. The aim of the UK method is to predict Z(x) at a nonsampled area as well. It splits the random function into a linear combination of deterministic functions, the smoothly varying and non-stationary trend, that is also called a drift m(x) e R, and a random component Y(x): = Z(x) À m(x), representing the residual random function [10]. OK assumes a stationary, i.e. constant mean of the underlying real-valued random function Z(x). But in reality, the mean value varies, it is not often constant across the entire study area and the variable seems to be non-stationary. A non-stationary regionalized variable can be considered with two components [11]; drift (average or expected value of the regionalized variable) and a residual (difference between the actual measurements and the drift). The method of UK assumes that the mean m(x) has a functional dependence on the spatial location and can be approximated by a model with the equation [12]: where a l is l th coefficient to be estimated from the data, f l is l th basic function of spatial coordinates that describes the drift, and k is the number of functions used in modeling the drift. The OK and UK functions used in this study are including: stable, hole effect, J-Bessel, and Gaussian.

Model evaluation
After the interpolation, each spatiotemporal point in the check data set would have both the original PM 2.5 and PM 10 measurements and an estimated value. Then we conducted four accuracy assessments to compare the original and estimated PM 2.5 and PM 10   O i s are the original values [13]. The accuracy assessments are defined as follows (Eqs. (10)-(13)): The MAPE, also known as Mean Absolute Percentage Deviation (MAPD), is a measure of prediction accuracy of a forecasting method in statistics, for example in trend estimation. It usually expresses accuracy as a percentage.
The RMSE has been used as a standard statistical metric to measure model performance in meteorology, air quality, and climate research studies.
The MBE statistical indicator is also commonly used in comparing the models of predictions. Low values of MBE are desirable, but overestimation of an individual data element will cancel underestimation in a separate observation.
The MAE is another useful measure widely used in model evaluations. The MAE is a quantity used to measure how close forecasts or predictions are to the eventual outcomes.

Statistical analysis
Statistical description of data (PM 2.5 and PM 10 ) for interpolation techniques was performed by XLSTAT. The spatial structure analyzing data and zooning by interpolation methods include: Ordinary Kriging, Universal Kriging, Inverse Distance Weighting, and Radial Basis Functions done by Arc GIS software version 10.3. Finally, the best interpolation method was chosen according to the values of each error algorithm. It should be noted that in case of normal data, interpolation methods will have the best results. In case of non-normal data distribution, we can do normalization by "Cox-box". The statistical analysis for geographical comparisons was done by the SPSS software (Chicago, IL, USA, Version 16.0). The data were tested for normality using a Shapiro-Wilk test. The data were not normally distributed. We used non-parametric procedures, and Kruskal-Wallis test followed by a Mann-Whitney U test. A multiple regressions was calculated between particulate matter and independent meteorological parameters. The Spearman correlation test was used to examine the relationship between PM 2.5 and PM 10 particulates as well as between particulate matter and meteorological parameters.

Particulate matter standards
Based on standard of Iran Department of Environment (DOE) in 2014 which considered the annual maximum concentration for PM 2.5 as 10 mg/m 3 , the values obtained for the period of nine months were greater than the standard limit; but based on Environmental Protection Agency (EPA) standard in 2012 and standard of DOE, which announced the maximum concentrations of 24-h for PM 2.5 as 35 mg/ m 3 , the results were greater than standard in 75% of cases (Fig. 2). Regarding the standard of DOE in 2014, which announced the annual maximum concentration for PM 10 as 20 mg/m 3 , the values obtained for the period of nine months were greater than the standard limit, but compared with DOE 24-h for PM 10 as 154 mg/m 3 , the values were less than the strategy values in 89.58% of the cases. Furthermore, in comparison with EPA 24-h for PM 10 as 150 mg/m 3 , the values were less than the strategy values in 87.5% of the cases. In similar studies, the great amount of aerosol particles has been confirmed in other cities of Iran [14]. The ratio of PM 2.5 /PM 10 showed variability from 0.22 (Mashhad Road station) to 0.82 (Shahid Beheshti Square). It indicates that coarse particles (greater than 2.5 mm) make up the majority of aerosol (Table 1). In similar studies by USEPA, the annual mean PM 2.5 /PM 10 ratios measured in urban and semi-rural areas were between 0.3 and 0.7 [15]. The average PM 2.5 /PM 10 ratio during the sampling period was 0.50 compared to the range from 0.15 to 0.25, reported by EPA [16]. The percentage ratio of the mean concentrations of PM 10 /PM 2.5 /PM 1.0 in Tehran was found to be approximately as 7: 2: 1 [17].

Monitoring
The mean concentrations of particulate matter in 48 sampling stations during the ninth month study period for present research are shown in Table 1  in the Ghuchan Road station was the occurrence of strong wind entered from the east Sabzevar that it has significant impact on the air quality of the city. This station has a marginal position and is devoid of vegetation. Towhid station also has a high concentration. The station is located in the north Sabzevar. The lack of vegetation and construction activities could be due to high concentrations. The pollution in Sarbedaran station is caused by high traffic of heavy vehicles in the Mashhad highway located in this region. Furthermore, the main reason for the pollution in Police station, Sabzevar-Tehran is industrial activities in this region. In the south of Beyhaq Oil station, agricultural activities take place that can be a source of particulate matter.
In urban agglomerations, there are different and various point and line sources of particles. Whereas industrial activities such as domestic heating with coal or oil are seen as point sources, emission of particles from motorized traffic occurs mainly along the roads to be realized a line source. Emissions by motorized vehicles do not only include exhaust particles, but also involve abrasion products from tyres, brakes, clutches, and the road's surface. Furthermore, particles are emitted by resuspension of previously deposited particles by vehicle-induced turbulence. Besides the local emissions, particle concentrations in cities are also influenced by advection due to particle transport from rural surroundings or long-range, often trans-boundary transport. Moreover, the local winds seem to contribute to pollution in Sabzevar. In a study in India, it was shown that presence of dust and great traffic affected the amount of particles by 42.6% and 3.36% [18]. In Srilanka, the average of PM 2.5 concentration varied from 18 to 83 mg/m 3 in outdoor [19]. In a study in 6 municipal zones of Chile, the amounts of particles were greater in the central regions of the cities [20]. As shown in Table 1 and Fig. 3, the concentration of PM 10 is much higher than the PM 2.5 . It seems that PM 10 as pollutant is responsible for the pollution of Sabzevar. According to research conducted in USA [21], PM measurement results showed that highway and marginal urban areas by improving traffic bottleneck, the amount of particulate matter is reduced to 41 percent. The studies of on the mass PM 2.5 and PM 1 effects in Helsinki urban air pollution showed that PM 2.5 is the most effective PM index for air pollution which is most significantly associated with respiratory and Cardiovascular disease [22]. Measurements of PM 10 and PM 2.5 in urban area of Nanjing, China showed that more than 70% of total suspended particles are of a size that they are deposited in the respiratory tract below trachea, whereas about 22% of the mass is respirable and will reach the alveoli [23].

Spatial modeling
In this study, the data from 48 stations were normalized by applying Box-Cox method and the coefficient of 0.6 for PM 2.5 and 0.1 for PM 10 . The Box-Cox transformation is a particularly useful family of transformations. Normality assumptions are critical for many univariate intervals and hypothesis testing. It is important to test the normality assumption. If the data are in fact not normal, the Box-Cox normality plot can often be used to find a transformation that will approximately normalize the data. The distributions and normality of the observed data and the predicted values using the RBF, IDW, OK, and UK techniques were analyzed by using the histograms, the Q-Q plots, and the box plots  Fig. 4, RBF histogram shows an example close to a normal distribution; OK and UK histograms show a distribution slightly right skewed; while IDW histogram have normal distribution (p-value of Kolmogorov-Smirnov test = 1.000). Fig. 5 shows that the p-value of Kolmogorov-Smirnov test for IDW is 0.997 and OK box plot shows a distribution highly skewed to the right. These results indicate that the performance of IDW is slightly better than the other methods.
It is worth mentioning that in the present study the mean concentration during the 9-month study period were used in order to draw the distribution map. After applying interpolation function based on the above-mentioned method, pollution zoning for each pollutant (PM 2.5 and PM 10 particulate  Fig. 8). Eastern winds are prevalent during study period. In the official meteorological reports, long term records also show that the prevailing wind direction has been east, owing to the topography of the region. Looking at the maps of these methods (Figs. 6 and 7), the pollution of particulate matters is observed in the eastern regions of Sabzevar, which may be due to the prevailing wind blowing from the east in Sabzevar. The pollution is also observed in the southwest caused by traffic could be the result of high traffic of heavy vehicles in the highway located in this region (Mashhad highway), which is the entrance to Sabzevar. Figs. 6 and 7 show the generated maps of the PM 2.5 and PM 10 predicted using RBF, IDW, OK and UK for the data sets in Sabzevar. The Figures show a difference between predicted maps using the deterministic (RBF and IDW) and geostatistical methods (OK and UK) techniques. Fig. 6a and b shows more red zones than c and d. Fig. 6b shows how the IDW algorithm behaves when a high PM outlier is next to a low PM value; the abrupt change in prediction surface is not seen with other methods. From Fig. 6b, the IDW outperformed in estimating values from the observed data. Fig. 7 is the same way. In this (Figs. 6 and  7), the contribution of sampling density is very significant. The large sampling density makes the performance of OK and UK worse than the performance of IDW and RBF. However, it is not enough to evaluate the model performance only through the estimated values and maps. Comparing the estimated values with the observed values was used to evaluate the performance of all IDW, RBF, OK and UK techniques in terms of the accuracy of estimates [24]. In this study, the comparison of performance between the interpolation techniques was achieved by using RMSE and MAPE. The lower value of MAPE and RMSE for each interpolation method indicates the optimality of that method. Generally, lower the estimated error, the interpolation method previously mentioned to prepare the [ ( F i g . _ 6 ) T D $ F I G ] maps of pollution zoning are more appropriate. Table 2 presents the comparison of MAPE and RMSE values for each method of PM 2.5 and PM 10 particulate matters. Moreover, two parameters of MBE and MAE were used. The MAE and the RMSE were used to evaluate the accuracy of predictions. The MAE is similar to the RMSE but less sensitive to large forecast errors. The RMSE is the square root of the variance of the residuals. It indicates the absolute fit of model to the data and how close the observed data points are to the model's predicted values [25]. Whereas R-squared is a relative measure of fit, RMSE is an absolute measure of fit. As the square root of a variance, RMSE can be interpreted as the standard deviation of the unexplained variance, and has the useful property of belonging to the same units as the response variable. Lower values of RMSE indicate a better fit. As Table 2 shows, the value of RSME for IDW method is 0.023 (PM 2.5 ) and 0.03733 (PM 10 ), respectively which is the lowest value among the methods. Examining MBE revealed some points. First, the model estimates the value of the variable more or less. Second, it determined the value quantity. When MBE is equal to zero, it indicates that the model has estimated the value under investigation in a good way without any diversion. The accuracy of the model is determined by MAE parameter. Regarding MAE, the value of zero is the indication of 100% accuracy of the model, and more distance from zero reflects lower accuracy of the model. This method of investigation was applied for all models of zoning. In general, the results shown in Table 2 indicate some differences in the values of MAPE, RMSE, MBE, and MAE. This means IDW was better than other methods in estimating the values from observed data. Karydas et al. [4] showed that the applied interpolation methods have similar results in terms of accuracy without clear advantages than to each other. Furthermore, in another study [26], no differences were found between Kriging and IDW methods for some elements, and MBA, MAE and RMSE factors were similar for these two methods; but for other studied elements, Kriging estimator showed a higher accuracy. Halek el al. [17] estimated urban suspended particulate in Tehran. The model was made based on the data from 42 stations located within the region. To build a "Surface Model" for PM 10 , PM 2.5 and PM 1.0 , different algorithms were used to interpolate the data from those obtained for the known site. Then, the results were extended to the "surface". For this purpose, mean concentrations of PM 10 , PM 2.5 and PM 1.0 , in each sampling site for June, July and August were calculated, interpolated and generalized to the surface by "Inverse Distance Weight" or "Spine" algorithms by using ArcGIS 9.2. The results strongly indicated that the concentrations of PM 10 , PM 2.5 and PM 1.0 of any points inside the region, including the traffic zone, fail to meet the required international standard values. The extracted estimate concentrations for the 22 hospitals reveal that the concentration of PM 10 for "Azadi Psychic", "Children" and "Mustafa Khomeini" hospitals are the worst, estimating from the model to be 119.42 mg/m 3 , 107.09 mg/m 3 and 101.14 mg/m 3 respectively.
Understanding spatial variability of air pollutant concentrations appear to be critical for public health assessments. Bermana et al. [27] examined ground-level ozone and evaluated the performance method for predicting and mapping national concentrations across the United States, while assessing the significance of accounting for spatial variability. Ozone concentration was predicted by using four approaches, including Land Use Regression (LUR), IDW, OK, and UK, and evaluated with a Monte Carlo cross-validation simulation. Results were mapped for the continental United States. UK was preferred over OK by allowing us to assess the significance of environmental covariates both for inference and prediction of ozone concentrations. In another research [28], two different methods, namely as Kriging method and Inverse Distance Weighted method, were examined for developing Digital Elevation Model image. Each method's advantages and disadvantages were considered. Here, the Kriging yielded better estimates.

Multiple regressions of particulate matters
A multiple regressions was calculated between particulate matter and independent meteorological parameters. Furthermore, the Spearman correlation test was used to examine the relationship between particulate matter and meteorological parameters. All the results indicated a significant relationship. The relationship of PM 2.5 and meteorological parameters was positive for temperature (p = 0.046 and Pearson coefficient = 0.25) and wind speed (p = 0.006 and Pearson coefficient = 0.60) and negative for relative humidity (p = 0.001 and Pearson coefficient = À0.223) and precipitation (p = 0.03 and Pearson coefficient = À0.10). Furthermore, the relationship of PM 10 and meteorological parameters was positive for temperature (p = 0.046 and Pearson coefficient = 0.23) and wind speed (p = 0.006 and Pearson coefficient=0.58) and negative for relative humidity (p = 0.001 and Pearson coefficient = À0.221) and precipitation (p = 0.03 and Pearson coefficient = À0.11). Considering the significant relationship between meteorological parameters and the pollutants, meteorological variables were found to be effective in the air pollution of Sabzevar. The results of this study revealed that the highest concentration of particulate matter has occurred at high temperatures. In other words, with increased temperature, the concentration has increased as well, since the correlation between particulates and temperature is positive and significant. According to the results of statistical surveys, the maximum concentrations of pollutants have occurred in the low relative humidity rate. In other words, increased air dryness has been associated with an increase in the amounts of pollutants. Increased relative humidity, if accompanied by a rainfall phenomenon, can reduce the air pollutants through washing. Therefore, this is one of the main factors in reducing the amounts of pollutants during high relative humidity and rainfall. There is a positive correlation between the mean concentration of particulate matter and the  [29]. The Spearman correlation test was used to examine the relationship between PM 2.5 and PM 10 particulates. The variations of 24 h simultaneous PM 10 and PM 2.5 concentrations were correlated with R = 0.76 (p = 0.000). This finding showed the strong correlation between PM 10 and PM 2.5 concentrations.

Additional information
Environmental monitoring can be defined as the systematic sampling of air, water, soil, and biota in order to observe and study the environment, as well as to derive knowledge from this process. Monitoring can be conducted for a number of purposes, including use of data to control the environment [30][31][32], to establish environmental "baselines, trends, and cumulative effects", to test environmental modeling processes, to educate the public about environmental conditions, to inform policy design and decisionmaking, to ensure compliance with environmental regulations, to assess the effects of anthropogenic influences, or to conduct an inventory of natural resources [33][34][35]. Air pollution appears to be one of major dilemmas of urbanization and crowded cities. Particulate matter from the perspective of public health seems to be a major air pollutant. Concerns over ambient fine particulate matter (PM) pollution are increasing in the modern world due to its potentially harmful effects on the human health and the environment [36]. Dust may lead to climate changes on a global scale and local, as well as inducing changes in biological cycles and the environment. For, air quality is a major factor in the public health, which depends on the concentrations of particulate matter. This has been confirmed by comparing PM concentrations with life expectancy [37]. PM 1.0 , PM 2.5 and PM 10 are particulate matter with an aerodynamic diameter less than 1.0, 2.5 and 10 mm, respectively. Particles smaller than 10 mm (PM 10 ), as a life expectancy index, include small liquid droplets and solid particles that can easily be inhaled deeply [38]. Once threshold PM concentration levels exceed over prolonged periods, particulate matter might can a wide range of disturbances and illnesses [39]. Inflammation of the eyes, lungs, and the skin are other adverse effects of personal PM exposure [40]. Due to harms to the health of people caused by exposures to air pollutants in urban areas, monitoring and forecasting of air quality parameters have become popular as an important topic in atmospheric and environmental research today [41].
In air pollution studies, the air quality models are used to predict and estimate the concentrations of one or more species within the space and time relevant to the dependent variables. Modeling enables us to assess the current as well as the future air quality to make "informed" policy decisions [42][43][44][45]. Interpolation methods have been well developed to estimate values at unknown locations based upon values spatially sampled in GIS. They are characterized as either deterministic or stochastic depending on using statistical properties [46]. Deterministic models include Inverse Distance Weighted (IDW), Rectangular, Natural Neighbours, and Spline methods. The IDW interpolator assumes that each input point has a local influence, which diminishes with distance. It gives more weight to the points closer to the processing cell than those are in further distant. Geostatistics is the second method used, consisting of Kriging in its various forms, Ordinary Kriging (OK), Simple Kriging (SK), and Universal Kriging (UK) [47,48]. Kriging is a geostatistical technique to estimate the values of random fields at unobserved points resulted from observing values at known locations. IDW and Kriging techniques are widely used interpolation techniques [4]. Kravchenko and Bullock [46] evaluated the effect of data variability and the strength of spatial correlation in the data on the performance of the grid soil sampling of different sampling density and two interpolation procedures, Ordinary Point Kriging and Optimal Inverse Distance Weighting. Wollenhaupt et al. [49] compared these two interpolation methods and concluded that the IDW is more accurate for zooning P and K levels of soil. Setianto and Triandini [28] assessed the Kriging interpolation method to be more accurate compared to the IDW. Berman et al. [27] evaluated many methods for spatial mapping and came to the conclusion that Kriging was more accurate for zoning O 3 levels of air ambient. Summarizing the above, the accuracy of each method depends on the assumptions and subjective judgments made, such as using or not-using the results smoothness and linearity of interpolation functions. In a study, the concentrations of PM 10 , PM 2.5 and PM 1.0 were measured in urban areas of Tehran at warm and cold seasons, which data were applied in the related modeling using the Arc-GIS. To this end, the samples were collected from 42 sites in an 18 km 2 region located in the west and central parts of Tehran. The mean concentrations of PM 1.0 , PM 2.5 and PM 10 were respectively reported as 13.14 mg/m 3 , 22.67 mg/m 3 and 95.72 mg/m 3 in the warm season; and 50.12 mg/m 3 , 70.72 mg/m 3 and 193.86 mg/m3 in the cold season [50]. Shad et al. [51] measured the particulate matter with a mass median aerodynamic diameter of less than 10 mm with concentrations at 52 sample stations in Tehran to identify dangerous areas for the human health. Followed by PM 10 data prediction, their results demonstrated that genetic algorithms can reduce the estimated error (3.74) compared to linear functions (8.94 and 12.29). In recent years, atmospheric models, such as GIS, have been used for environmental analysis and the related management for supporting the environmental decision makers in different countries. In this study, concentrations of PM 2.5 and PM 10 were found in urban areas of Sabzevar for the first time and the data were applied in the related modeling by using Arc-GIS. This study was an attempt to evaluate spatial interpolation methods to determine the concentration of PM 2.5 and PM 10 in Sabzevar and select the most suitable interpolation method for preparation of zoning maps of the particulate matter in GIS. Multiple regressions analysis were used to investigate the effect of independent variables (meteorological parameters) on particulate matters concentrations.
In this research, the mean concentrations for PM 2.5 and PM 10 were as 56.61 mg/m 3 and 110.93 mg/ m 3 , respectively. According to the results of this study, the high pollution of particulate matters is observed in the eastern regions of Sabzevar which may be due to the prevailing wind blowing from the east in Sabzevar. The high pollution is also observed in the southwest caused by traffic could be the result of high traffic of heavy vehicles in this highway, as an entrance to Sabzevar. The results from this study showed that concentrations of PM 2.5 and PM 10 are higher than the regulatory standards in most stations; this calls for urgent need for continuous monitoring and control of source point emission as well as proper awareness for the citizenry for health hazards associated with air pollution. The data distribution and the autocorrelation were set as parameters to evaluate RBF, IDW, OK, and UK techniques. The best interpolation method for the particulate matter (PM 2.5 and PM 10 ) was deterministic method by IDW function. In summary, one can say that the model is useful for the estimation of pollutant concentrations within urban areas.