Statistical Character and Transport Pathways of Atmospheric Aerosols in Belgrade

Clean air is considered to be a basic requirement for human health and well being. Various chemicals are emitted into the air from both, natural and anthropogenic sources. In spite of the introduction of cleaner technologies in industry, energy production and transport, air pollution remains a major health risk and tighter emission controls are being enforced by many governments. Atmospheric particles – aerosols – are some of the key components of the atmosphere. They influence the energy balance of the Earth's surface, visibility, climate and environment as a whole [1-3]. According to World Health Organization (WHO), ozone, particulate matter (PM), heavy metals and some hydrocarbons present the priority pollutants in the troposphere [4]. Public health can also be indirectly affected by deposition of air pollutants in environmental media and uptake by plants and animals, what results in entering of chemicals into the food chain or drinking water, and thereby constituting additional sources of human exposure. A number of epidemiological studies have demonstrated that acute and chronic health effects are related to the inhalable PM10 (aerodynamic diameter less than 10 μm) exposure in the urban environment, and some data also seem to indicate possible seasonal effects of the particulate matter on human health [5-10]. This is especially important for urban aerosols, whose variety of size and composition make complete characterization a difficult task. Particulate matter pollution is nowadays one of the problems of the most concern in great cities, not only because of the adverse health effects, but also of the reducing atmospheric visibility and affect to the state of conservation of various cultural heritages [11]. Therefore, the measurement of the levels of atmospheric particulate matter is a key parameter in air quality monitoring throughout the world.


Introduction
Clean air is considered to be a basic requirement for human health and well being. Various chemicals are emitted into the air from both, natural and anthropogenic sources. In spite of the introduction of cleaner technologies in industry, energy production and transport, air pollution remains a major health risk and tighter emission controls are being enforced by many governments. Atmospheric particles -aerosols -are some of the key components of the atmosphere. They influence the energy balance of the Earth's surface, visibility, climate and environment as a whole [1][2][3]. According to World Health Organization (WHO), ozone, particulate matter (PM), heavy metals and some hydrocarbons present the priority pollutants in the troposphere [4]. Public health can also be indirectly affected by deposition of air pollutants in environmental media and uptake by plants and animals, what results in entering of chemicals into the food chain or drinking water, and thereby constituting additional sources of human exposure. A number of epidemiological studies have demonstrated that acute and chronic health effects are related to the inhalable PM10 (aerodynamic diameter less than 10 μm) exposure in the urban environment, and some data also seem to indicate possible seasonal effects of the particulate matter on human health [5][6][7][8][9][10]. This is especially important for urban aerosols, whose variety of size and composition make complete characterization a difficult task. Particulate matter pollution is nowadays one of the problems of the most concern in great cities, not only because of the adverse health effects, but also of the reducing atmospheric visibility and affect to the state of conservation of various cultural heritages [11]. Therefore, the measurement of the levels of atmospheric particulate matter is a key parameter in air quality monitoring throughout the world.
As a result of health and environmental impacts of PM, more rigorous regulations are in force in the USA and European countries. PM standards, issued by European Commission (EC), have included PM10 monitoring and limit values in the Air Quality Directive in 1999 [12]. Directive established in the first stage, annual limit value of 40 μg m -3 and 24h limit value of 50 μg m -3 (not to be exceeded more than 35 times in a calendar year) to be met by 2005, and in the second stage annual limit value of 20 μg m -3 and 24h limit value of 50 μg m -3 (not to be exceeded more than 7 times in a calendar year) to be met by 2010. The detailed discussion of these limit values, regulations and relations of new EU standards to US EPA standards can be found elsewhere [13,14].
Most of the trace metals are emitted in particulate form [15][16][17] and present in almost all aerosol size fractions, but are mainly accumulated in the smaller particles [18]. This has a great effect on the toxicity of metals since the degree of respiratory penetration depends on particle size. Trace metals are persistent and widely dispersed in the environment and interacting with different natural components result in toxic effects on the biosphere [19,20]. The studies of the transport and mobilization of trace metals up to now have attracted attention of many researchers [21][22][23][24]. Trace metals are released into the atmosphere by human activities, such as combustion of fossil fuels and wood, high temperature industrial activities and waste incinerations. The combustion of fossil fuels constitutes the principal anthropogenic source for Be, V, Co, Ni, Se, Mo, Sn, Sb, and Hg. It also contributes to anthropogenic release of Cr, Mn, Cu, Zn, and As. High percentages of Ni, Cu, Zn, As, and Cd are emitted from industrial metallurgical processes. Exhaust emissions from gasoline may contain variable quantities of Ni, Cu, Zn, Cd, and Pb [25,26]. In addition to the PM mass concentration limit values, also based on health impact criteria, recent European Union standards set target (Ni, As, Cd) and limit (Pb) values for metals [12]. Environmental technologies may have to be adopted in specific industrial spots to reach the target values. For aimed reduction of PM10 levels and their trace metal contents detailed knowledge of sources and their respective contribution to the PM10 levels, is required.
One of the main difficulties in air pollution management is to determine the quantitative relationship between ambient air quality and pollutant sources. Various receptor models have been used to identify aerosol sources and estimate their contributions to PM10 and trace elements concentrations at receptor sites [27][28][29][30][31][32]. The combined application of different types of receptor models could possibly solve the limitations of the individual models, by constructing a more robust solution based on their strengths. Processes in the atmosphere represent a complex problem due to the simultaneous influence of several independent factors. Similarly to other air pollutants, PM10 concentrations are random variables influenced by the emission level, meteorological conditions and topography. Each area is specific and the required emission reduction to meet Air Quality Standard (AQS) is different. Information about the frequency distribution of pollutants is useful for developing air pollution control strategies. When the specific probability function of an air pollutant is known, it is easy to predict the required emission reduction, the frequency of exceedance of the AQS, as well as the return period [33]. Many types of probability distributions have been used to fit air pollutant concentrations. These include lognormal distribution [34,35], Weibull distribution [36], gamma distribution [37] and type V Pearson distribution [38]. The Weibull probability model is suitable for presenting air quality data in a "pseudolognormal" distribution [39,40]. In the present study, daily average mass concentrations of PM10 were taken from 2003 to 2005 in order to estimate the parameters for three theoretical distributions. Lognormal, Weibull and type V Pearson distributions were chosen to fit the PM10 data in Belgrade. The distributional parameters were estimated by the maximum likelihood method. Based on the fitted distributions the minimum reduction required in order to achieve AQS was estimated using rollback equation [41][42][43]. Moreover, since the tail of theoretical distributions diverged in the high concentration region, a two-parameter exponential distribution was used to fit high PM10 concentrations and estimate exceedances and return period more precisely [44]. These methods can reasonably predict the return period and exceedances in the succeeding period and can also be used for predicting source emission reduction.
In many countries, the weekly variation in the concentration of atmospheric aerosols and meteorological variables in urban areas has been reported as evidence of anthropogenic influence on the weather system. Of the primary sources road transport is the dominant contributor to high particulate level in the air in most major cities. This study examines the temporal and weekly variations in meteorological variables and PM10 mass concentrations in Belgrade. The cumulative frequencies and quantile-quantile plot for Sundays against workdays were constructed. The results showed that Sunday concentrations were significantly reduced from those during the week.
The Mediterranean region and particularly the Balkan Peninsula have been under the influence of Saharan dust transport and deposition over millennia. Identification of the concentration, composition, origin, transport and geographical distribution of PM10 in Mediterranean atmosphere has been the subject of research activities since the last two decades as it is heavily affected by two contrasting sources; mineral dust (mainly from Sahara Desert) [45] and various anthropogenic (from industrialized/semi-industrialized countries) emissions. Deviation from lognormal distribution of atmospheric aerosols is a positive indicator of possible transport processes [33]. Since we are particularly interested in transportation of trace metals content in PM10 samples, the quantiles for some trace elements were analyzed by the quantile-quantile P-P slope test. This was used for fitting the quantiles of the expected theoretical cumulative probability for normal, lognormal and Weibull distribution with the quantiles of the cumulative probabilities calculated for the experimental dataset [46].
In addition, we investigated the transport pathways and potential sources of PM10 concentrations based on backward trajectories and PM10 concentration records from 2004 to 2006. Cluster analysis was used to reveal the major pathways for different seasons as well as corresponding statistical analysis related to different clusters. Hybrid receptor models Potential Source Contribution Function (PSCF) and Concentration Weighted Trajectory (CWT) were used for identification of source regions. The PSCF values can be interpreted as a conditional probability describing the spatial distribution of probable geographical source locations by analysis of trajectories arriving at the sampling site [47]. Since the PSCF method is known to have difficulties in distinguishing strong from moderate sources, the CWT model that determines the relative significance of potential sources has been additionally performed [29,48].
In this review, we report some of the results of the integral monitoring of air quality in Belgrade urban area in order to evaluate the impact of airborne trace metals on the pollution load for the period from 2003 to 2006. Some of the results concerning suspended particle mass and trace metal concentrations in ambient air of the Belgrade will be presented, with the aim to examine elemental associations and to indicate the main sources of trace and other metals in the city. The results of this long-term project of the pollution monitoring could be used as the baseline data for analysis of health risks due to inhalation of suspended aerosols, and to provide scientific evidence for setting up an air pollution control strategy. This information is crucial in environmental quality assessment, and can lead to the determination of a possible exceedance of the critical loads. This point can be considered as traffic-exposed.

Experimental methods and procedures
The second sampling point was on the roof of the Rector's Office building of Belgrade University, at a height of about 20 m. As this sampling point is in the very city centre, on the rooftop where the airflow is not blocked by any direction, it can be considered as representative for urban-background concentrations. Suspended particles were collected on pre-conditioned and pre-weighed Pure Teflon filters (Whatman, 47 mm diameter, 2 μm pore size) and Teflon-coated Quartz filters using MiniVol air sampler (Airmetrics Co. Inc., 5 l min -1 flow rate) provided with PM10 cutoff inlet. Particulate matter mass concentration was determined by weighting of the filters using a semi-micro balance (Sartorius, R 160P), with a minimum resolution of 0.01 mg. Loaded and unloaded filters (stored in Petri dishes) were weighed after 48 hours conditioning in a desiccator, in the clean room at a relative humidity of (45-55)% and a temperature of (20 ± 2) 0 C. Quality assurance was provided by simultaneous measurements of a set of three ''weigh blank'' filters that were interspersed within the pre-and post-weighing sessions of each set of sample filters and the mean change in "weigh blank" filter mass between weighing sessions was used to correct the sample filter mass changes. After completion of gravimetric analysis, PM samples were digested in 0.1 N HNO3 on an ultrasonic bath. An extraction procedure with dilute acid was used for the evaluation of elements which can become labile depending on the acidity of the environment [49]. This procedure gives valid information on the extractability of elements, since the soluble components in an aerosol are normally dissolved by contact with water or acidic solution in the actual environment.
The elemental composition (Al, V, Cr, Mn, Fe, Ni, Cu, Zn, Cd, and Pb) of the aerosol samples was measured by the atomic absorption spectroscopy method (AAS). Depending on concentration levels, samples were analyzed for a set of elements by flame (FAAS) (Perkin Elmer AA 200) and graphite furnace atomic absorption spectrometry (GFAAS) using the transversely-heated graphite atomizer (THGA; Perkin Elmer AA 600) with Zeeman-effect background correction. The THGA provided a uniform temperature distribution over the entire tube length, rapid heating and an integrated L'vov platform, which gave an improved signal/interference ratio and high analytical sensitivity. Analyte injection (20 μl) and the atomization were done in five steps controlled by the appropriate software and auto-sampler. For calibration, standard solutions containing all metals of interest were prepared using Merck certified atomic absorption stock standard solutions containing 1000 mg l -1 metal in 0.5 N HNO3 and Milli-Q quality deionised water, with no matrix modifier addition. Details on sampling procedures and PM analysis are given in detail elsewhere [50][51][52][53][54][55].

Estimation of emission source reduction
Assuming unchanged spatial distribution of emission sources, meteorological conditions and no reactive species, according to the rollback equation [36], the emission source reduction R (%) required to meet AQS can be estimated as: In this equation E{C}s is the mean (expected) concentration of the distribution when the extreme value equals Cs (the concentration of the AQS), E{Cp} is the mean concentration of the actual distribution and Cb the background concentration. If Cs = 125 μg m -3 , the PM10 daily average concentration is not exceeded more than once per year (P[PM10>Cs] = 1/365 = 0.00274), then E{C}s is the expected daily PM10 average concentration of a distribution where the probability of a concentration exceeding 125 μg m -3 is equal to 0.00274.
Frequency distribution of PM10 is a very useful tool to estimate how frequently a critical concentration level is exceeded. In order to obtain the best frequency distribution of PM10, three theoretical functions were used: lognormal, Weibull and type V Pearson distributions. The distribution parameters were estimated by the maximum likelihood method. Probability density functions as well as corresponding equations for estimating distribution parameters are presented in Table 1.

Distribution
Probability distribution function Maximum likelihood estimated parameter The appropriateness of each distribution was assessed by the Kolmogorov-Smirnov (K-S) test. The K-S statistic is defined as the maximum difference between the sample cumulative probability and the expected cumulative probability i.e.
where ( ) n f x and ( ) F x are the expected and observed cumulative frequency functions, respectively. The D value will be compared to the largest theoretical difference, Dα/2 acceptable for the K-S test at a certain significance level α. If D is less than Dα/2, one can conclude that the hypothesis that there is no difference between observed and expected distributions has been accepted at the α level of significance. Generally speaking, if D value declines the goodness of fit increases.
The fitted results of three theoretic distributions and the measured data for PM10 are presented in Figure 1. Once the distribution parameters were determined, K-S test could be used to assess which type of distribution is more appropriate for representing the PM10 distribution. From Figure 1 it is clear that the Weibull distribution is inappropriate for representing PM10 distribution (D = 0.073), while type V Pearson distribution is the most suitable one (D = 0.038).   statistical textbooks. It was found that the value of E{C}s (the average concentration of lognormal distribution where the probability of a concentration exceeding 125 μg m -3 equals 0.00274) was 32.1 μg m -3 . Therefore, the mean PM10 concentration should be reduced from the observed value of 68.4 μg m -3 to 32.1 μg m -3 and the estimated emission reduction can be calculated from equation 1 being equal to 53.1%. For type V Pearson distribution the estimated source reduction required to meet AQS is 58.6%. Therefore, the source emissions should be controlled much more to reduce PM10 concentrations and meet AQS in the future period [41].
Previous theoretical distributions can give good result for estimating the mean concentration and required reduction of source emission. However, the fitted results of the parent distributions are not accurate enough in the high concentrations region. Therefore, a two-parameter exponential distribution was applied to predict the return period and exceedances of the critical PM10 concentration.
A two-parameter exponential distribution derived from extreme value theory [44] represents the cumulative frequency distribution of high concentrations over a specific percentile.
where , where N1 is the size of the chosen high PM10 concentration and 1 rN P is the probability of a value that is ranked r out of N1 values. The relation between variate n y and 1 rN P is: and the parameters bm and  can be estimated by the least-squares method. In addition, the return period ( ) c R x , defined as the average number of averaging periods (or observations) between exceedances of a given critical concentration, xc , can be calculated as: where f is the chosen specific percentile. After determining return period the days of exceedances within one year can be calculated in order to develop air control strategy.

Estimating the frequency distribution of PM10 by the statistics of wind speed
Since the distributions of air pollutant concentrations are influenced by meteorological conditions, such as wind speed, Simpson et al. [59] linked air quality data and wind speed to the Atmospheric Turbulence and Diffusion Laboratory (ATDL) model. The concentration of air pollutant data at cumulative probability, p, is inversely proportional to the wind speeds v at probability (100-p) when the frequency distributions of both data are log-normally distributed and the geometric standard deviation (shape factor) is the same e.i.
where K depends on the average area source strength of pollutant and atmospheric stability factor (which varies with atmospheric conditions). This model was used to predict the maximum level of air pollutant. Finally, the frequency distribution of air pollutants can be estimated from the frequency distribution of wind speed [60]. Here we demonstrate the relationship between the frequency distributions of wind speed and daily PM10 measured in Belgrade during 2003 year. Three theoretical frequency distributions, (lognormal, Weibull and gamma) were used to fit the measured PM10 data. The distributional parameters were estimated by the method of maximum-likelihood and the K-S test was used to determine which type of distribution is appropriate to represent the frequency distribution of PM10 and wind speed. The fitted results of PM10 are shown in Figure 4 and it can be seen that the lognormal distribution is the most appropriate one to represent the statistical character of PM10 (Dmax = 0.06). The fitted results of lognormal distributions for PM10 and wind speed, with corresponding parameters, are presented in Figure 5. If is found that the shape factors of lognormal distributions σg are similar, therefore equation 8 can be used to estimate the distribution of PM10 from the distribution of wind speed. The average K value can be calculated from the regression line using method of least-squares for the wind speed data and PM10 concentrations from the 30-th to 80-th percentiles (the K value is almost constant in that range). The obtain regression is  When the atmospheric condition is unchanged and the K value is estimated, the concentration of PM10 at different percentiles can be obtained. Observed and estimated frequency distribution of PM10 is shown in Figure 6. The difference between the observed and estimated concentration is significant in the higher cumulative probability range, because the inconsistency of the both, wind speed and PM10 data with lognormal distribution in this region. The result of K-S test shows that the difference between observed and estimated distribution is not significant at the 5% significance level, so the obtain results are acceptable. It is clear that there is no universal distribution that can be used to represent the frequency distribution of PM10 concentrations. This analysis shows that if the distribution type and shape factors are the same for the wind speed data and air pollutant concentration, there exists simple relationship between two data sets and the frequency distribution of the air pollutant can be reasonably estimated. It can be useful for understanding statistical characteristics of air quality in different region.

Weekend effect and PM10 temporal variations
Weekly cycles of the mass concentrations of anthropogenic aerosols have been observed in many regions and the intensity of these cycles vary significantly depending on region and season [61,62]. In addition, ozone concentrations tend to be higher on weekends than on weekdays, despite the fact that lower emissions of ozone precursors are expected on weekends. This phenomenon is known as the "weekend effect" [63]. The mechanisms for weekend effect on ozone formation are still not well understood. In a lot of papers that are dealing with weekend effect in detail, the diminution of NOx, CO, PM10 and volatile organic compound (VOCs) concentrations can be noticed during the weekend [64,65]. Although the relationship between the variation in aerosol concentration and weather variables has not been clearly understood, many studies have revealed significant weekday-weekend differences in weather variables such as the diurnal temperature range [66], and clouds and precipitation [67,68].
Generally speaking, several source types make significant contribution to particle level in Belgrade: combustion processes, mixes road dust, traffic, metal processing and transport [30,50]. Like in many other major cities with a few millions inhabitants, of the primary sources, fossil fuel combustion and traffic contributions play the most significant role [51].
Since traffic and commercial pattern in urban area show a marked change at weekends, it is expected that particulate concentrations will decrease. Markovic et al. [69] reported the mean daily concentrations of SO2, NO2, PM10, CO and O3 during the weekdays and weekend in November/December 2005 case study in Belgrade urban area (Figure 7). The results explicitly show that daily average concentrations of SO2, NO2, PM10 and CO show significant weekly cycles with the largest values around midweek and smallest values in weekend.
With the exception of ozone, concentrations of airborne pollutants were higher during weekdays. It has been suggested that the weekend condition of lower particle concentrations resulted in less absorption of sunlight, which in turn leaded to enhanced photochemical reaction and enhanced ozone formation [70]. Besides, the higher concentrations of PM10, NOx, CO and SO2 during the weekdays influence on the increased deposition of O3 on the particles as well as on the increased consummation of O3 in the reactions with NOx, CO and SO2 what also can have an impact.
In order to investigate differences in PM10 mass concentrations between Sundays and weekdays quantile-quantile (Q-Q) plots were constructed [71]. Q-Q plots are useful in detection the distribution differences in data set and allow comparison of concentration distribution of the entire range of values, from the smallest to the largest. If two cumulative frequency diagrams (x and y) are plotted, then, corresponding to any frequency p, there are two quantile values and and Q-Q plot is a scatter plot of versus for various p. Obviously, if x and y are identically distributed, the plot of versus would fall on a straight line with the slope equal to one.   Figure 10. It can be seen that the maximum concentrations occur in the morning and late evening which can be significantly related to traffic pattern and anthropogenic influence. Significant seasonal variations of some metals (Ni, V, Al and Fe) have already been reported [30] indicating the impact of fossil fuel combustion and local heating units. Monthly variation of PM10 mass concentrations for the year 2005 is presented in Figure 11 showing significant differences during the winter and summer period. Analysis of air quality data in the frequency domain contributes to the understanding of periodic behaviours and yields information about temporal and spatial scales of the underlying mechanisms [72]. The results obtained in many studies show that the frequency analysis of air quality time series is a powerful technique that can provide important information about the nature of the processes behind the measurements. This methodology could be applied to determine representative background concentrations of air pollutants by removing short-term fluctuations associated with influence of local emission sources [73].

Transport analysis
Processes in the atmosphere represent a complex problem due to the simultaneous influence of several independent factors (meteorological conditions, pollutant emission level and topography) and thus the environmental data are random variables that follow natural lognormal distribution [33]. Deviation from lognormal distribution of atmospheric aerosols is a positive indicator of possible transport processes. In this study the quantiles for some trace metals in PM10 samples were analyzed by the quantile-quantile P-P slope test. It is used for fitting the quantiles of the expected theoretical cumulative probability for normal, lognormal and Weibull distribution with the quantiles of the cumulative probabilities calculated for the experimental dataset.
Fittings of quantiles of the experimental database with the expected quantiles for theoretical normal, lognormal and Weibull distributions for Al, Mn and V are shown in Figure 12.
The P-P slope test results show that the Weibull distribution is the most appropriate for V and Ni suggesting their possible regional transport. The other metals followed lognormal distribution (except Cu and Zn data which are not successfully described by any of theoretical distribution function) what is a positive signal that an active sources exist in the vicinity. For unknown metal sources, hybrid models that incorporate wind trajectories Potential Source Contribution Function (PSCF) and Concentration Weighted Trajectory (CWT) were used to resolve source locations [29].

Potential Source Contribution Function (PSCF)
Air parcel back trajectories, ending at the receptor site, are represented by segment endpoints. Each endpoint has two coordinates (latitude, longitude) representing the central location of an air parcel at a particulate time. To calculate PSCF, the whole geographic region of interest is divided into an array of grid cells whose size is dependent on the geographical scale of the problem so that PSCF will be a function of locations as defined by the cell indices i and j. If a trajectory end point lies at a cell of address (i, j), the trajectory is assumed to collect material emitted in the cell. Once aerosol is incorporated into the air parcel, it can be transported along the trajectory to the receptor site. The PSCF value can be understood as a conditional probability that describes the spatial distribution of probable source locations. The staying time of all trajectories in a single grid cell is nij, and mij is the staying time in the same cell that corresponds to the trajectories that arrived at a receptor site with pollutant concentrations higher than a pre-specified criterion value (average PM10 and metal concentrations were used in this study). The PSCF value for the ij-th cell is then defined as Cells related to the high values of potential source contribution function are the potential source areas. However, the potential source contribution function maps do not provide an emission inventory of a pollutant but rather show those source areas which emissions can be transported to the measurement site. To remove the large uncertainty caused when a grid cell has small staying time and large PSCF values, usually a weight function W(nij) multiplies the PSCF value to better reflect the uncertainty in the values for these cells [27]. In this study the grid covers area of interest defined by (35 0 -50 0 ) N and (10 0 -30 0 ) E with cells 0.5 0 x0.5 0 latitude and longitude.
Air masses back trajectories were computed by the HYSPLIT (HYbrid Single Particle Lagrangian Integrated Trajectory) model [74] through interactive READY system [75].
The weighting function reduced the PSCF values when the total number of the endpoints in a particular cell was less than about two times the average value of the end points per each cell.

Concentration Weighted Trajectory (CWT)
In the current PSCF method, grid cells having the same PSCF values can result from samples of slightly higher concentrations than defined by the criterion or extremely high concentrations. As a result, larger sources can not be distinguished from moderate sources. Therefore, a method of weighting trajectories with associated concentrations (CWTconcentration weighted trajectory) was developed. In this procedure, each grid cell gets a weighted concentration obtained by averaging sample concentrations that have associated trajectories that crossed that grid cell as follows: Similar to the PSCF method, CWT method also employs the arbitrary weight function to eliminate grid cells with few endpoints. Weighted concentration fields show concentration gradients across potential sources. This method helps to determine the relative significance of potential sources.
For air mass trajectory visualization and statistical analysis, a software application called TrajStat was used in which clustering, PSCF and CWT methods were included [76].

PSCF and CWT Results
Some studies have already reported a long-range transport of PM from western countries over Balkan Peninsula which is sporadically (mostly in spring and summer) associated with African dust outbreaks [77,78]. Based on the trace metals data set contents in PM10 for the period 2003-2005 the possible transport process over Belgrade was investigated. The PSCF and CWT plots for PM10, V, Al and Mn are presented in Figure 13. It can be seen that PSCF plot for PM10 shows high probability area in the northeast and west region. CWT plot indicated that the high concentrations of PM10 are related to the sources located in the west, northwest and northeast direction distinguishing major sources from moderate ones by calculating concentration gradients. The highest PSCF and CWT values for V were similarly distributed in the northeast region. Since V followed Weibull distribution, transportation process is expected and this region the most likely influences receptor site. Aluminum and Mn are metals related to dominant local sources. As can be seen from Figure 13, there are no significant emission region that Al can be transported from, while Mn transportation process can be expected from southeast direction.

Clustering of modelled backward trajectories
During 2004-2008, 48-hour backward trajectories of air masses over Belgrade were daily calculated at the previously defined six different altitudes. Transport paths of air masses with similar history and origin were obtained by cluster analysis. Since we were concerned on the directions of the trajectories, the angle distance between back trajectories [76] has been used as the cluster model. The angle distance between two backward trajectories was defined as The variables x0 and y0 define the position of the studied site and 1 1 ( ), ( ) x i y i and 2 2 ( ), ( ) x i y i are coordinates of i-th segment for trajectories 1 and 2, respectively. In order to investigate possible seasonal variation in transportation process of PM10 over Belgrade, we analysed summer and winter season separately. Cluster analysis resolved the four main directions of air mass flows over Belgrade during both periods. A cluster number can be deduced through visual inspection and comparison of the mean-trajectory maps. The representative trajectories for each cluster group are presented in Figure 14. As can be seen, the four classes were related to west southwest, east, west northwest and north directions. Trajectories within different classes of airflow had distinct effects on the PM10 concentrations. The mean PM10 concentrations and frequency of air back trajectories occurrences within each cluster group are presented in Figure  15. During the summer period the most frequently arriving direction of air masses is from west northwest direction (40%) and evenly distributed among the others cluster groups. Despite the high frequency of occurrences, the contribution from this direction to PM10 concentration measured at the receptor site is the least. High PM10 concentrations are associated with the trajectories coming from west southwest and east directions. Also, the highest number of air back trajectories from these directions is associated with the observed PM10 concentrations above critical level of 50 μg m -3 . During the winter period the most frequently arriving direction of air masses is from west southwest direction (37%) and north direction (26%). Although the frequency of air back trajectories coming from east is the lowest, their contributions to the high PM10 concentrations are the most significant. The highest average contribution to the observed PM10 is from the cluster group representing the arrival direction from west southwest. In these directions there are several possible emission sources: coal-fired thermal power plants -"Nikola Tesla A and B", Obrenovac (west southwest); "Kostolac A and B", Kostolac (east southeast); "Kolubara", Veliki Crljeni (south southwest); as well the steel industry complex in Smederevo (east southeast).
Obviously this kind of analysis can reveal the different transport pattern of particulate matter during summer and winter season. Cluster analysis of air back trajectories can be useful tool for investigation of the major transport pathways and their contributions to the high PM10 values. Additional analysis can be used to estimate background and long distance transportation of particulate matter [79,80].

Conclusion
The atmospheric processes are very complex, and pollutant concentrations are inherently random variables because of their dependence on the fluctuation of a variety of meteorological conditions and emission intensity. When sets of air pollution data are available and we aim to predict future trends, we have to know some characteristics of the measured variables, therefore, various statistical characteristics can be determined (maximum and minimum concentration, standard deviation, frequencies of exceedances of some critical level, etc) and assigned to the pollutant concentration. The distribution function is one of the most important characteristic of the random variables. Although there is a priori no reason to expect that "atmospheric" distribution should adhere to a specific probability distribution in some cases, there may be a natural probability distribution which relates emission levels to pollution levels. A great number of papers dealing to this problem have been published and in this study three theoretical distributions (lognormal, Weibull and type V Pearson) were used to fit PM10 concentrations in the Belgrade urban area while two-parameter exponential distribution was used to fit the high percentile PM10 daily concentration region. Based on the fitted distributions the minimum reduction required to attain the daily AQS was estimated using the rollback equation. These methods can reasonably predict the return period and exceedances in the succeeding period. The results suggested that the required PM source reduction ranged from 53% to 63% in Belgrade area. Recently measurements indicated that some progress in this direction has already been done and the PM10 mass concentrations have been reduced during the last few years.
In addition, it is demonstrated that the distribution of PM10 can be successfully estimated from the wind speed data if the corresponding two distributions followed the same type of distribution and have similar shape factors.
Daily and monthly variations of PM10 were analyzed and confirmed the existing seasonal variation. Weekend effect has also been observed and the differences between weekdays and weekend were analyzed by Q-Q plot and cumulative distribution function.
The analysis of quantiles of real and theoretical distributions for trace metals in PM10 indicated dominant local emission sources for the most of metals. To identify the source locations associated with the high fluxes of trace metal contents in PM10 (Al, Mn and V), two trajectory-based models, PSCF and CWT, were used. Both PSCF and CWT resolved similar source locations for those elements. The results indicated that the source areas contributing to receptor site were located mostly in north east and west directions.
In addition, trajectory clustering was further used to investigate the transport pathways during summer and winter season. Four clusters were generated for both periods and corresponding average contributions to PM10 mass concentrations. Although the frequency of air back trajectories coming from east was the lowest during the summer period, their contributions to the average PM10 concentration were the most significant (41.7 μg m -3 ) followed by the trajectories coming from west south west (40.4 μg m -3 ). The highest average contribution to the observed PM10 concentration during the winter season was from the cluster group representing the arrival direction from west southwest (74.4 μg m -3 ). Obviously this kind of analysis can reveal the different transport pattern of particulate matter during summer and winter season.
Nowadays commercial instruments with online measurement are helpful and can provide a large base of measurement. Still, researchers have to analyses these measurements and mathematical models and statistical analysis are some of the useful methods for developing air pollution control strategy related to specific geographic region.