ARTICLE IN PRESS
Atmospheric Environment 39 (2005) 3661–3674 www.elsevier.com/locate/atmosenv
Comparative assessment of regionalisation methods of monitored atmospheric deposition loads Frido Reinstorfa,, Maja Binderb, Mario Schirmera, Jost Grimm-Strelec, Wolfgang Waltherd a
UFZ–Centre for Environmental Research Leipzig-Halle Ltd., Department Hydrogeology, Theodor-Lieser-Str.4, D-06120 Halle, Germany b RWTH Aachen University, FIW, Mies-van-der-Rohe Str. 17, D-52056 Aachen, Germany c State Office for Environmental Protection Baden-Wu¨rttemberg, PSF 210752, D-76157 Karlsruhe, Germany d Dresden University of Technology, Institute for Groundwater Management, Karcherallee 7, D-01277 Dresden, Germany Received 9 August 2004; received in revised form 28 February 2005; accepted 4 March 2005
Abstract The objective of this investigation is to assess the suitability of well-known regionalisation methods of data from existing deposition monitoring networks for use in water resources management. For this purpose a comparison of the applicability and accuracy of various regionalisation methods was made. A crucial point is the data demand of the various methods. In this investigation the deterministic and geostatistical methods inverse distance weighting (IDW), ordinary kriging (OK) and external drift kriging (EDK) as well as the chemical transport models METRAS-MUSCAT, EMEP, EDACS and EUTREND have been characterised and evaluated. The methods IDW and OK have been applied to the investigation areas—the German Federal States of Lower Saxony and Saxony. An evaluation of these methods was carried out with a cross-validation procedure. The result was in most cases a higher accuracy for the OK method. The EDK method has been investigated in order to find suitable drift variables from the parameters precipitation amount, altitude and wind direction. With help of a correlation analysis a suitable drift variable could not be found. After the application of OK, verification was carried out by a comparison of the estimated data set with an independently determined data set. The result was a relatively smaller deviation of the estimated data set. The investigation considers data from routine monitoring networks as well as networks for special applications and has been carried out on the basis of monitoring networks of the two states. The investigated database was wet and bulk + 2+ 2 deposition of the substances NH+ , and Cd2+ in Lower Saxony and SO2 4 , SO4 , NO3 , Na , Pb 4 in Saxony. From this, a consistent database of bulk deposition data was built. From all applied methods OK proved to cope best with the data deficiencies that were found. r 2005 Elsevier Ltd. All rights reserved. Keywords: Deposition loads; Regionalisation; Geostatistics; Monitoring networks; Deposition models
1. Introduction
Corresponding author. Fax: +49 345 5585 559.
E-mail address:
[email protected] (F. Reinstorf).
With regard to the protection and sustainable development of soils and waters, there is a growing need for the authorities to obtain integrated and area-specific
1352-2310/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.atmosenv.2005.03.006
ARTICLE IN PRESS 3662
F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
environmental information. Typically, the measurement of atmospheric deposition and substance transport in soils and groundwater is carried out on the basis of sectoral monitoring networks, i.e. networks for deposition, seepage water and groundwater. Because each network was designed without regard to the others, the monitoring stations of the various networks are located independently. For this reason methods have to be developed, that can be used to link available information at various monitoring points in various monitoring networks. Common methods for the regionalisation of seepage and groundwater data are geostatistical and deterministic methods. The objective of this investigation was to find a comparable method to regionalise available monitoring data from existing routine deposition networks with a maximum accuracy and acceptable effort for water resources management. The project was initiated to exemplify a suitable regionalisation method for the networks in the German Federal States of Lower Saxony and Saxony. The concept was to assess well-known regionalisation methods as well as chemical transport models (CTMs) with regard to their accuracy and feasibility under given conditions. Special networks from various institutions should be considered. The investigated database was for a time period from 1995 to 1998. The investigation was carried out by means of bulk deposition data.
2. Monitoring networks, database and data treatment In Germany, an official countrywide deposition monitoring network is operated by the Umweltbundesamt (UBA). It consists of 37 measuring stations (status quo in 2001) for the monitoring of regional pollutant transport. Furthermore, the German Federal States operate monitoring networks to fulfil different research objectives, e.g. air hygiene, documentation of input paths or registration of basic impacts of emissions far from the source. Additional operators of networks on local and regional scales exist (see Gauger et al., 1997, 2000, 2002). The number of available monitoring stations amounts to 76 in Lower Saxony (area: 47.343 km2; i.e. one station per 623 km) and 24 in Saxony (area: 18.338 km2; i.e. one station per 764 km2). The locations of the stations are shown in Fig. 1. The observed substances are, in general, NH+ 4 , NO3 , + + 2+ 2+ + 2+ 2+ 2 SO4 , Cl , Na , K , Mg , Ca , H , Cd , Pb and the amount of precipitation. For this investigation 2 the following substances were selected: NH+ 4 , SO4 , + 2+ 2+ NO3 , Na , Pb and Cd in Lower Saxony and SO2 4 in Saxony. The selection was carried out on the basis of the needs formulated by the water resources stakeholder group of the project and the quantity and quality of the
available database. The total number of monitoring stations that could be used for the investigation period 1995 to 1998 were as follows. In Lower Saxony, 44 stations could be used for the substances ammonium, sulphate and sodium, 57 for nitrate and 18 for cadmium and lead. In Saxony, 19 stations could be used for sulphate. The sampling techniques in Germany that are applied in the Federal States may be divided into wet-only and bulk-sampling. Therefore, wet-only measurements have to be converted to bulk. This was done according to Gauger et al. (2000). Thereby, the wet-only were increased by a parameter-dependent factor (f p , see Table 1) based on the following equation: X p-Bulk ¼ X pWet 1=f p ,
(1)
where X pBulk is the calculated amount of the parameter p in the bulk sampler, X pWet is the measured amount of the parameter p in the wet sampler and f p is the factor according to Table 1 for parameter p [–]. The parameter p is the acronym of the substance to be measured. Sulphate data were corrected with the so-called ‘‘sea-salt correction’’ after Gauger et al. (1997). Because rainwater contains ions originating from seawater that is transported into the continent by sea spray the sea-salt correction is used to correct the measured ion content from a sample for the sea-salt contribution. This sea-salt influence would produce a trend in the data set which overlays the impacts of other sources. Eq. (2) is based on the assumptions that sodium is 100% sea-borne near the coast and the ratio of element concentrations in seawater and sea spray is equal. These assumptions are plausible for northern and western Europe (Gauger et al., 1997) as well as for the wet deposition in East Germany (UN ECE, 1996). SSSC ¼ Sdep Nadep
S SW , NaSW
(2)
where S ssc is the corrected sulphur deposition in consideration of sea spray emissions (SSC ¼ sea-salt corrected) in [moleq ha1 a1], Sdep is the deposition of sulphur in [moleq ha1 a1], Nadep is the deposition of sodium in [moleq ha1 a1] and SSW/NaSW is the concentration ratio of sulphur to sodium in sea water which is equal to 0.120 (UN ECE, 1996). The data of the networks are available in different time steps (weeks, monthly averages or sums). At some monitoring stations the time series are incomplete. Especially in winter there are some months without any data, e.g. because of plugging of the samplers by snow. In order to guarantee a consistent database, the missing data were interpolated. For this, in some cases
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3663
Fig. 1. Deposition observation networks in Lower Saxony and Saxony and general map of Germany.
up to 6% of the monthly values had to be interpolated. With this, a database of the yearly deposition sum was built by multiplying the mean concentration of a year with the yearly precipitation amount. A monthly calculation was not performed because often consecutive months were without data. Therefore, we considered
that the data series were too uncertain to use for the comparison of methods. Applying Eq. (3) (Gauger et al., 2002), the quality of the measured field deposition is checked by calculating the ion balance from the main components in a liquid water phase. The assumption is that the amounts of the
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3664
Table 1 Results of the ratio wet to bulk-deposition (f p ) from simultaneous measurements (after Gauger et al., 2002, supplemented with data from Van Leeuwen et al., 1996)
Mean (f p ) Standard deviation
NH+ 4
SO2 4
NO 3
Na+
Pb2+
Cd2+
0.89 0.17
0.89 0.17
0.83 0.11
0.86 0.20
0.71 0.20
0.73 0.14
0.85 0.09
0.85 0.04
0.86 0.07
— —
— —
For comparison: Van Leeuwen et al. (1996) Mean 0.83 Standard deviation 0.1
main cations and anions in solution have to be in equilibrium:
IonBalance
¼
et al., 1996) is a ‘‘3D, non-hydrostatic meteorology model’’ that works with a ‘‘ground following co-
2þ ðNHþ þ Mg2þ þ Kþ þ Naþ þ Hþ Þ ðSO2 4 þ NO3 þ Cl Þ 4 þ Ca 100, 2þ 2 þ 2þ þ þ þ ðNH4 þ Ca þ Mg þ K þ Na þ H Þ þ ðSO4 þ NO3 þ Cl Þ
where IonBalance is the ion balance in (%) and the ions have the unit [eq]. If ion balance errors exceeded 20%, the data sets—that means all data (see Eq. (3)) measured at one sampling campaign at one sampler—were rejected. Less than 2% of the data sets had to be rejected. For the application of the methods, there needed to be at least 30 monitoring stations with consistent data sets available. This could be fulfilled for all substances except for the heavy metals Pb2+ and Cd2+ in Lower Saxony (only 18 stations) and for SO2 4 in Saxony (only 19 stations).
3. Investigated regionalisation methods The investigated methods are the deterministic method of inverse distance weighting (IDW), and the geostatistical methods of ordinary kriging (OK) and external drift kriging (EDK). Furthermore, processbased chemical transport model (CTMs) were evaluated. A brief description of the methods follows. 3.1. Chemical transport models Different models for various aims and spatial scales are available to regionalise substance depositions from the atmosphere. Three models for typical applications were selected and are explained briefly in the following. The coupled model system METRAS-MUSCAT was built for the simulation of reactive substances in the meso-scale and was applied in the INTERREG—II— Project OMKAS (OMKAS, 2000) in order to estimate the effects of selected emission reducing measures in the so-called ‘‘Black triangle’’ (the boundary region of Germany, Poland and the Czech Republic) on the immissions and depositions. METRAS (Schlu¨nzen
(3)
ordinate system’’. MUSCAT (Knoth and Wolke, 1998) is a CTM, which is able to simulate the transport, deposition and transformation of substances by using the chemistry mechanism Euro-RADM (Stockwell and Kley, 1994). Considered boundary conditions are, among others, the topography and the land use in high spatial resolution (OMKAS, 2000). The long-range transport model EMEP (European Monitoring and Evaluation Programme for Investigations of the transborder transport of air pollution) (Erisman and Draaijers, 1995) belongs to the deterministic Lagrange models and is used for the quantification of the loads of the substances SO2, NH3, NO, NO2, 2 HNO3, NO and NH+ 3 , SO4 4 , which are transported transborder. In order to determine a total deposition map (150 150 km raster), the immission concentrations are calculated from emission maps (50 50 km grid). The model EDACS (European Deposition of Acidifying Components on a Small Scale) estimates the dry deposition for six land-use classes in Germany using the inferential method (van Leeuwen et al., 1996; Draaijers et al., 1996). Within EDACS the EUTREND model (van Jaarsveld, 1995; Bleeker et al., 2000) is used for assessing receptor-specific ammonia/ammonium (NHx) dry deposition fluxes. The substances SO2 and SO2 4 , aerosols (SOx), NO, NO2 (NOx), HNO3 and NO 3 + aerosols, NH3 and NH+ 4 aerosols (NHx) and Na , Ca2+, K+, Mg2+ aerosols are considered. The dry deposition of the base cations Na+, Mg2+, Ca2+ and K+ is determined from regional annual wet deposition (results of mapping the wet deposition in the 1 1 km grid) with the help of scavenging ratios. The latter are computed with a single leaching model which is based on simultaneous measurement of cations in the surfacenear air and in precipitation (mean of annual concentrations of a number of measurement points). The
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
EDACS mapping output is in 1 1 km grid resolution for the land-use classes urban areas, water surfaces, agricultural areas, coniferous, deciduous and mixed forests derived from CORINE (Coordination of Information on the Environment) land cover map. Finally, the total deposition can be estimated in a grid of 10 20 km, by interpolation of measured data of wet deposition overlaid by dry deposition, which is calculated using a resistance model, e.g. EDACS. 3.2. Deterministic interpolation methods The basis of deterministic interpolation methods is the principle of simple averaging. The accuracy of the results of these methods is, solely influenced by the available data density and distribution. The IDW method weights the data to estimate a non-observed point, ui, inversely proportional to the squared distance di from the estimation point (Isaaks and Srivastava, 1989): , Xn 1 Xn 1 zðu Þ , (4) zIDW ðuÞ ¼ i 2 i¼1 d i¼1 d 2 i i zIDW ðuÞ
is the estimated value at point u, zðui Þ is where the observed value of variable z at monitoring point ui , d i is the distance of monitoring point i to the point to be estimated, and n is the number of monitoring points. In contrast to simple arithmetical averaging this method considers the distance between points. An increasing distance between a monitoring point and a point of estimation means a decreasing influence of the monitoring point. 3.3. Geostatistical interpolation methods A large number of different methods belong to the group of kriging methods, which are included in geostatistical interpolation methods. As a common characteristic they consider the spatial distribution of measurement values by using variogram analysis. For a stationary random function a variogram is defined as follows: 2gðhÞ ¼ VarfZðu þ hÞ ZðuÞg, gðhÞ ¼ Cð0Þ CðhÞ; 8u,
ð5Þ
where CðhÞ is the stationary covariance, and Cð0Þ ¼ VarfZðuÞg is the stationary variance. More generally the variogram is a measure of the spatial variability and the variogram distance measures the average degree of dissimilarity between an unsampled value zðuÞ and a nearby data value. Ordinary kriging (OK) (e.g. Matheron, 1971, cited in Bardossy, 1993; Deutsch and Journel, 1998) is a univariate method, estimating a value at a point of a region for which a variogram is known. Thereby data of the neighbourhood are used. OK is an exact interpolator
3665
because at the monitoring points the estimated and observed values are identical. In other words, the simulated parameter distributions are conditioned. Stationarity of the random function is assumed. The OK-estimator is zOK ðuÞ ¼
n X
lðOKÞ ðuÞzðui Þ, i
(6)
i¼1
ðuÞ where zOK ðuÞ is the estimated value at point u, lðOKÞ i are the weights and zðui Þ is the observed value of variable z at monitoring point ui : More accurate results than with univariate OK can usually be obtained by multivariate geostatistical methods such as EDK (Ahmed and De Marsily, 1987). EDK allows the consideration of transient conditions. The EDK-estimator is zEDK ðuÞ ¼
n X
lðEDKÞ ðuÞzðui Þ; i
and
(7)
i¼1 n X
lðEDKÞ ðuÞC R ðuj ui Þ þ m0 ðuÞ þ m1 ðuÞyðui Þ j
j¼1
¼ C R ðu ui Þ; a ¼ 1; :::; n, n X lðEDKÞ ðuÞ ¼ 1; j j¼1 n X
lðEDKÞ ðuÞyðuj Þ ¼ yðuÞ; j
ð8Þ
j¼1
’s are the kriging (EDK) weights, the m’s where the lðEDKÞ i are Lagrange parameters and C R is the covariance. However, a well-correlated drift variable must be available for the main variable. In this way, the information density can be enhanced. In EDK, the spatial-variate external drift variable Y must be a linear function of the spatial-variate expectation E of the main variable z at point u. The trend model is then EfZðuÞjY ðuÞg ¼ a þ b Y ðuÞ,
(9)
where E is the expected value, ZðuÞ is the parameter value of the additional information at location u, a and b are parameters and Y ðuÞ is the spatial-variate external drift variable.
4. Assessment of the regionalisation methods 4.1. Applicability of chemical transport models The advantage of CTMs is their ability to simulate deposition with high resolution in time and space. For example, EDACS output is a grid with a relatively high resolution, compared with deposition monitoring networks and a more or less high degree of process-relation.
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3666
In contrast to that, the disadvantage is the extensive amount of input data required, e.g. EDACS uses 6-h averages of the meteorological parameters radiation, temperature, air moisture, wind speed and others. Furthermore, EDACS is not based on observed deposition data, but on immission concentrations simulated by the EMEP model. The model EMEP—a large-scale module—uses the chemical composition of precipitation on a daily basis at 106 observation stations in 25 European countries. The model METRAS uses emissions in the best resolution available for the species to be calculated. In addition, initial values for all the different species are needed. The implemented chemical reaction scheme is, however, either the compact EMEP scheme with 11 reactions and 13 species, or a photolytic scheme with 173 reactions and 73 species. This means the same time resolution mentioned above with regard to EMEP is needed. These facts show there are manifold reasons why the models described are not feasible for regionalisation of available deposition data from monitoring networks in Germany: 1. The deposition data from existing monitoring networks are currently not available in the required time resolution. 2. The necessary parameters are hard to obtain. 3. The data processing requires significant effort.
4.2. Comparison of interpolation methods For the application of OK and EDK variogram modelling has to be performed (see Section 4.2.1). The OK and IDW tests were followed by EDK (see Section 4.2.2). The software used was STATGRAPHICS for data analysis, ArcView (ESRI) for the creation of maps and for IDW and GEO-EAS (Englund and Sparks, 1991) for OK. 4.2.1. Variogram modelling For kriging methods optimal variogram models were fitted on the basis of calculated experimental variograms. The fit is optimal in the sense that the weighted residual sums of squares are minimised. The variogram model was chosen from the following list: exponential, spherical, linear, circular and Gaussian. In most cases the spherical model was optimal (Table 2). Fig. 2 presents the experimental variograms with the models fitted. 4.2.2. Criteria for the comparison of methods For comparison the cross-validation procedure is suitable. It permits an objective comparison of the methods. To assess the variogram model, after Davis (1987), cited in Sinowski (1995), the standard deviation between the interpolated values and the observed values at the monitoring points is compared with the calculated mean kriging standard deviation. Sðui Þ ¼
Furthermore, these models are very sophisticated, meaning only a few scientists are able to use these models. In addition, it should be mentioned that CTMs are originally not designed for regionalisation of depositions but for process-based understanding and simulation. However, if the total deposition is the aim there is, in our view, no alternative for the use of CTMs to calculate the dry deposition, which has to be added to the interpolated wet deposition.
standard deviation , s0 ðui Þ
(10)
where standard deviation is the standard deviation between the interpolated and observed values and s0 ðui Þ is the calculated mean kriging standard deviation. A variogram model can be accepted, if Sðui Þ is between 0.9 and 1.1. This criterion is suitable for data sets originating from micro-scale investigations, e.g. soil characteristics at an investigation plot size of 10 10 m. For depositional data this is hard to meet. Therefore, the best fit
Table 2 Variogram models of the depositions (4-year-average: 1995–1998) and the parameters Substance
Model
Nugget
Sill
Range (km)
Variance
Lower Saxony SO2 4 NO 3 Na+ NH+ 4 Pb2+ Cd2+
Exponential Spherical Spherical Spherical Spherical Spherical
0.1 3.2 0.4 1.3 135.0 0.02
13.9 10.8 8.3 3.5 410.0 0.08
161 60 145 150 155 80
13.8 7.6 7.9 2.2 275.0 0.06
Saxony SO2 4
Spherical
12
76
130
64
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
Sulphate 5
12 8 4
4 3 2 1 0
0 0
50000
(1)
Nitrate
50000 10000 150000 200000 250000 300000 Distance Sodium
16 12
12
Variogram
Variogram
0
10000 150000 200000 250000 300000 (2) Distance
16
8 4
8 4 0
0 0
40000
(3)
0
80000 120000 160000 200000 Distance (4) Cadmium
0.16
40000
80000 120000 160000 200000 240000 Distance Lead
600 500
0.12
Variogram
Variogram
Ammonium
6
Variogram
Variogram
16
3667
0.08 0.04
400 300 200 100
0 (5)
0
0
40000 80000 120000 160000 200000 240000 Distance (6)
0
50000
10000 150000 Distance
200000
250000
Sulphate 90
Variogram
75 60 45 30 15 0 (7)
0
50000
100000 Distance
150000
200000
Fig. 2. Variograms for deposition loads for the Federal States Lower Saxony (1) Sulphate, (2) Ammonium, (3) Nitrate, (4) Sodium, (5) Cadmium, (6) Lead and for Saxony (7) Sulphate. The distance is given in [m].
that could be obtained was accepted. Large variances may be the result of ion concentrations varying over shorter distances than those between monitoring points. To resolve this source of spatial variability, a denser and more accurate network is required.
Cross-validation can also be used to assess the various interpolation methods. In this case actual data are dropped one at a time and re-estimated from some of the remaining neighbouring data. The cross-validation should be as difficult as the actual estimation of
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3668
unsampled values. For example, if data are aligned along drill holes or lines, then all data originating from the same hole or line should be ignored when reestimating a datum value in order to approach the sampling density available in actual estimation (Deutsch and Journel, 1998). Each datum is replaced in the data set once it has been re-estimated. A continuous new estimation of the variograms was not carried out because this leads to the problem identified by Isaaks and Srivastava (1989), in that the simultaneous adjustment of all of the model parameters, including the type of function used, has not yet been proven successful in practice. The mean differences between observed and estimated values (Bias) and the mean squared prediction error (MSPE) are calculated. The main criteria for the suitability of the interpolation method are:
no systematic error (Bias is zero or close to zero), errors are minimal (MSPE as small as possible), and the coefficient of determination (COD) (Ha¨rdle and Simar, 2003) should be maximal (close to 1).
The equations used for the error calculations are the following: Bias ¼
n 1X ðz ðui Þ zðui ÞÞ, n i¼1
Table 3 Comparison of the errors for IDW and OK
(11)
MSPE ¼
n 1X ðzðui Þ z ðui ÞÞ2 , n i¼1
(12)
MSPE , (13) s20 where s20 is the total variance. In addition to the amount of emission and substance characteristics, the local exposition conditions may have an important influence on deposition. That could be additional information (drift variable) for EDK. The precondition is a correlation between the drift variable and the main variable. Possible influences on the following drift variables have been investigated: altitude, precipitation amount and wind direction (WD). In order to check the potential drift variables, a standard correlation analysis was carried out.
COD ¼
4.2.3. Results Table 3 shows the errors of the two methods OK and IDW for all investigated substances for Lower Saxony and Saxony. Sulphate is shown as yearly averages from 1995 to 1998 for Lower Saxony and 4-year-averages for Saxony. With regard to the cross-validation procedure for OK and IDW it can be shown for almost all substances, that Bias and MSPE, respectively, are smaller and COD is larger in most cases for OK (Table 3). COD permits the comparison of the results for various substances. For sodium achieved with OK, the maximum COD is 0.58.
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674 Table 4 Correlation matrix of the additional parameter altitude (ALT) and precipitation (P) for substance deposition in Lower Saxony and Saxony + SO2 NO Parameter ALT NH+ 4 4 3 Na
Lower Saxony ALT 1 0.07 0.55 P 0.59 0.34 0.75 Saxony ALT P
1 0.86
— —
0.68 0.84
0.57 0.27 0.81 0.03 — —
— —
Pb2+ Cd2+ P
0.07 0.73 — —
0.38 0.59 0.97 1 — —
0.86 1
The CODs of the other substances are much smaller, except for the heavy metals. As an example, sulphate also shows that for yearly averages OK provides smaller errors than IDW. In a nutshell, OK provides in most cases more accurate results than IDW. The standard correlation analysis applied on the additional parameters for EDK showed that altitude (ALT) is not correlated with deposition (main variable) for Lower Saxony (see Table 4). In contrast, precipitation amount (P) seemed to be relevant in Lower Saxony and in Saxony. But P is already used for the calculation of the substance deposition loads. Thus, the correlation is only ‘‘apparent’’ and P is not really additional information. The investigation of the dependency between wind director (WD) and deposition load was carried out on the basis of the wind fields, i.e. most frequent WD of all 38 monitoring stations in Lower Saxony. A relatively uniform wind field was found in Lower Saxony. This means the distribution of the most frequent WD does not vary significantly at different stations. Therefore, seven stations that are representative for the state area were selected for a more detailed investigation. The correlation analyses between WD and deposition load showed no correlations. This means different depositions can occur at one WD. Thus, this investigation also showed no significant correlation with regard to the deposition. Therefore, WD cannot be used as a drift variable for EDK. It can be concluded that useful additional information for the EDK procedure could not be found. But it should be mentioned here that WD could probably provide useful information in process-based atmospheric deposition simulations. In view of the comparison between OK, IDW and EDK it can be concluded that OK copes best with the available data. In the following, the verification of the data regionalisation with OK is explained.
3669
5. Verification Strictly speaking, verification of the OK data is not possible. This can only be done if the ‘‘true’’ values could be known. The only reliable information that can be used is the database of the monitored wet or bulk deposition loads. Because OK is a conditioned method, i.e. the OK data at the measurement points are always equal to the measurement data, a verification can only be made with help of another database. To accomplish this task, the database from Gauger et al. (2002) which was built by the UBA, can be used. Because this database does not contain bulk deposition data only wet deposition can be used for the verification. So, a minimum of about 70% of the total deposition is considered. The database will be denoted in the following as ‘‘UBA database method’’. The UBA database is split into three parts containing (a) precipitation monitoring, i.e. deposition data of samples analysed from wet and bulk samplers exposed in open field locations, (b) deposition data of bulk (wetonly and dry) throughfall and stemflow from deciduous forests, and (c) deposition data of bulk throughfall from coniferous forests. Bulk deposition, throughfall and stemflow fluxes are basic input for canopy budget model (CBM) calculations. The application of CBM allows calculation of estimates of dry, cloud and fog and total deposition for a forest stand by taking into account canopy exchange processes. An overview of canopy model calculations can be found in Draaijers et al. (1998). Wet and bulk deposition data, corrected to wet deposition fluxes, are used as input for mapping wet deposition fields. Throughfall and stemflow data are used to apply CBM calculations to derive dry and total deposition estimates, serving as evaluation data for modelled dry and total deposition maps. Verification was carried out with the wet deposition + + fields of SO2 for Lower Saxony 4 , NO3 , NH4 and Na 2 and SO4 for Saxony. It proceeded as follows. A regular grid of 50 50 km without consideration of the monitoring stations was constructed for both Federal States. For Lower Saxony the grid consists of 22 nodes and for Saxony 15 nodes. For every node a deposition was readout from the database. For these data, the basic statistical values mean value, standard deviation and percentage difference between OK and UBA method were calculated (Table 5). Furthermore, a correlation analysis was carried out (Fig. 3). Table 5 shows that both methods calculate similar mean values. The standard deviation has also equal magnitude. The percentage difference between the methods is smaller than 10%, except for Na+ with about 28%. In view of the relatively small differences one can say that the OK method generally produces a consistent database and the data estimated by OK are plausible in relation to the data of the UBA method.
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3670
Table 5 Comparison of the UBA method and the OK method: mean [kg ha1 a1] respectively for Pb2+ and Cd2+ [g ha1 a1], standard deviation and percentage difference [%] of the wet deposition for Lower Saxony and Saxony Substance
Lower Saxony SO2 4 NO 3 NH+ 4 Na+ Pb2+(b) Cd2+(b) Saxony SO2 4
Method
Mean of wet deposition
Mean of bulk deposition
Standard deviation
UBA OK UBA OK UBA OK UBA OK OK OK
23.21 21.55 21.20 19.51 9.39 10.14 10.63 13.57 31.65 1.36
26.1 24.1 25.5 23.4 10.5 11.4 12.4 15.7 44.6 1.9
4.73 3.7 3.32 3.02 1.91 1.49 10.04 10.39 12.94 0.51
UBA OK
19.06 20.95
21.4 23.5
3.26 4.73
Percentage difference to UBA(a)
7.2 8.0 8.0 27.7 — —
10.0
a Percentage difference to UBA ¼ (((Mean of wet deposition of the OK-method) (100))/(Mean of wet deposition of the UBA method))100. b UBA data are not available for the substances Pb2+ and Cd2+.
6. Discussion The suitability of the discussed regionalisation methods in water resource management depends always on the specific aim and the available database. In general, the issue of ‘‘suitable’’ or ‘‘not suitable’’ can only be addressed by discussing the advantages and disadvantages of specific methods. The suitability of specific methods for different objectives in water resource management is evaluated in the following. It can be ascertained that the monitoring networks of the German Federal States have been optimised with respect to an even distribution of locations within the state areas (Fig. 1), the sampling method and the kind of observed substances. This implies that the estimation of representative state-related deposition data can be achieved. In the case of the estimation of representative data for a sub-area of a state, the problem of data shortage may arise. With regard to the interpolation methods, first the deterministic procedure IDW should be discussed and after this OK and EDK. IDW enables the estimation of substance depositions from the atmosphere for non-monitored sites with small computational effort (UBA, 1995). The interpolation is independent of the scale of observation. The accuracy is influenced solely by the density and the distribution of the monitoring points. IDW is an interpolation technique in which interpolated estimates are made, based on values of nearby locations weighted only by distance from the interpolation location. IDW makes no assumptions about spatial relationships except for the
basic assumption that values at nearby points are more closely related to the interpolated location than values at distant points. Considering the density of monitoring networks, the amount of data and the high variability of the deposition in time and space, an application for smaller areas than Federal States and a higher time resolution than yearly averages cannot be recommended. With these restrictions, regionalisation on the basis of existing networks is possible. The accuracy is relatively high and acceptable. The geostatistical procedures, OK and EDK, have the advantage to consider the spatial correlation of the evaluated parameter. This is represented by a variogram. As demonstrated by the cross-validation analysis, more precise results can be achieved with the same database by a kriging procedure than by the abovementioned IDW method, provided that a spatial correlation of the considered parameter exists. This can be assumed, e.g. for the data with a high proportion of wet deposition, because at greater distances from the source the concentration decreases as a rule. It was demonstrated that the investigated database fulfils these conditions. Additionally, the spatial correlation consisting of the considered deposition parameter influences the accuracy of the interpolation results. In consideration of the available database which is theoretically too small to carry out a rigorous spatial analysis, e.g. variogram analysis, the application of OK on areas smaller than Federal States and on time steps smaller than one year
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
Nitrate (wet)
y=0.932x R2=0.3512 OK
OK
Sulphate (wet) 40 35 30 25 20 15 10 5 0
0
5
20 15 25 UBA-method
10
(1)
30
35
35 30 25 20 15 10 5 0
40
y=0.9328x R2=0.3615
0
5
10
30
35
y=0.8905x R2=0.7573
40 OK
OK
50
y=0.9834x R2=0.5225
25
Sodium (wet)
60
30 20 10 0
0
2
4
6 8 10 UBA-method
OK
(3)
15 20 UBA-method
(2)
Ammonium (wet) 16 14 12 10 8 6 4 2 0
3671
14
0
16
10
20
(4)
30 40 UBA-method
50
60
Sulphate (wet)
35 30 25 20 15 10 5 0
y=1.143x R2=0.5314
0 (5)
12
5
10
15 20 25 UBA-method
30
35
Fig. 3. Correlation between random grid wet depositions after UBA method and OK for Lower Saxony and Saxony: (1) Sulphate (wet) Lower Saxony, (2) Nitrate (wet) Lower Saxony, (3) Ammonium (wet) Lower Saxony, (4) Sodium (wet) Lower Saxony, (9) Sulphate (wet) Saxony.
cannot be recommended. For the nationwide determination of the wet deposition in Germany, the OK method has already been applied to estimate Critical Loads and Critical Levels (Gauger et al., 1997; Ko¨ble et al., 1997). More precise results than those obtained by the univariate OK method can be achieved by a multivariate geostatistical procedure, such as EDK, provided that a variable is available which correlates well with the considered parameter and that the variable can also increase the information density. Successful applications of this method can be found in references, e.g. for the interpolation of the amount of precipitation (e.g. Bardossy, 1993; Haberlandt and Kite, 1998), of soil water content (e.g. Lehmann, 1995) and of temperatures (e.g. Hudson and Wackernagel, 1994). With respect to the regionalisation of deposition data with geostatistical
procedures including additional information, no references were found. Until now CTMs have been used for more differentiated area-related estimations with the inclusion of additional information. In connection with substance deposition from the atmosphere, there are CTMs and/or link-ups available for different spatial and time scales. The advantage of a very high degree of detail description—both in time and space—is usually connected to the disadvantages of the amount of input data required, a low accuracy (Van Leeuwen et al. in Gauger et al., 1997) and the large amount of computational time needed. It should also be noted that these approaches are highly sophisticated and require the expertise of a specialist. Furthermore, the models, e.g. EDACS, do not use the monitoring data of the deposition networks;
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
3672
Emden
Emden Lueneburg Oldenburg
Hanover Braunschweig Osnabrueck
suiphate loads Yera average 1995-98 sea salt corrected
< 14 14.1 - 17.0 17.1 - 20.0 20.1 - 23.0 23.1 - 26.0 26.1 - 28.0 29.1 - 32.0 32.1 - 35.0 > 35
Goettingen
Osnabrueck
Ammonium Loads Year average 1995-98
Emden
Lueneburg
Nitrate Loads Year average 1995-98
Nitrate Loads [kg/ha*a]
Goettingen
Bremen
Hanover Braunschweig
Osnabrueck
< 15 15.1 - 18.0 18.1 - 21.0 21.1 - 24.0 24.1 - 27.0 27.1 - 30.0 < 30.0
Sodium Loads Year average 1995-98
Goettingen
(4)
(3)
Emden
Ammonium Loads [kg/ha*a] < 7.1 7.1 - 8.0 8.1 - 9.0 9.1 - 10.0 10.1 - 11.0 11.1 - 12.0 12.1 - 13.0 > 13.0
Lueneburg Oldenburg
Bremen
Hanover Braunschweig
Goettingen
(2)
Emden
Osnabrueck
Bremen
Hanover Braunschweig
Sulphate Loads [kg/ha*a]
(1)
Oldenburg
Lueneburg Oldenburg
Bremen
Sodium Loads [kg/ha*a] 3.4 - 10.5 10.6 - 17.6 17.7 - 24.7 24.8 - 31.7 31.8 - 38.8 38.9 - 45.9 46.0 - 55.0 > 55.0
Emden
Lueneburg Oldenburg
Lueneburg
Oldenburg Bremen
Bremen
Hanover Braunschweig Osnabrueck
codmium Loads Year average 1995-98
Hanover Braunschweig
Codmium Loads [g/ha*a]
Osnabrueck
0.7 - 0.9 1.0 - 1.2 1.3 - 1.5 1.6 - 1.8 1.9 - 2.1 2.2 - 2.4
Goettingen
(5)
Lead Loads Year average 1995-98
Goettingen
(6)
Lead Loads [kg/ha*a] 12.1 - 20.1 20.2 - 28.2 28.7 - 37.0 36.5 - 44.5 44.6 - 52.6 52.7 - 60.7 > 60.8
Leipzig Goerlitz Dresden chemnitz Sulphate Loads [kg/ha*a]
Zwickau
Sulphate Loads Year average 1995-98 sea salt corrected
< 15.0 15.1 - 19.0 19.1 - 23.0 23.1 - 27.0 27.1 - 31.0 31.1 - 35.0 > 35.0
(7) Fig. 4. Interpolation maps (OK) of the yearly averages of deposition loads [kg ha1 a1] for the Federal States Lower Saxony (1) Sulphate, (2) Ammonium, (3) Nitrate, (4) Sodium, (5) Cadmium, (6) Lead and for Saxony (7) Sulphate.
ARTICLE IN PRESS F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
rather immission maps are extracted from a simulation by the EMEP model. From these maps—partly with the help of other model approaches—the total deposition is determined. Other models, e.g. METRAS-MUSCAT, use monitoring data but implement relatively complex chemical reaction schemes which are hard to parameterise and the input data are mostly not available. It should be mentioned that if total deposition is the aim dry deposition models cannot be ignored. By comparing the different methods and data requirements, the OK method can be chosen as most suitable for regionalisation with acceptable accuracy and feasibility for water resource management issues. A set of maps of depositions for Lower Saxony and Saxony interpolated on the basis of this method (Fig. 4) have been created. The maps are made on a regular grid of 10 20 km. The means of wet and bulk deposition are shown in Table 5. As shown in the created maps (Fig. 4) the nitrate and ammonium deposition in Lower Saxony have similar emission patterns. These substance depositions, which are mainly caused by agricultural land use and animal farming, are correlated with the intensity of agricultural activity. While the highest nitrate deposition occurs outside the towns, the ammonium deposition is more homogeneously distributed. This is caused by the higher atmospheric mobility of ammonium gas. The high density of animal farms in the west of Lower Saxony can be shown very clearly by ammonium deposition. The nitrogen input of ammonium is about two times higher than nitrate. The sulphate deposition correlates also with the nitrogen deposition. It can be assumed that sulphate is mainly caused by agricultural sources, e.g. sulphatecontaining fertilisers, and by combustion of fossil burning. Furthermore, the deposition patterns of all mentioned substances show lower deposition near the coast of the North Sea. Sodium shows increasing deposition in a NW direction, caused by sea spray. Cadmium and lead depositions have ‘‘hot spots’’ in the northern part of Lower Saxony, correlating with the locations of harbours and shipyards. Also higher depositions are observed in the southern part. This correlates with the locations of the metal processing industries. The situation in Saxony is characterised by increasing sulphate deposition in the SE direction. The influence of emissions on the area coming from the Czech Republic is evident. Here, in the so-called ‘‘black triangle’’, situations are known where the sulphate concentration in the atmosphere is extremely high, especially in the winter months.
7. Conclusions The comparison of various regionalisation methods showed that the OK method is the best method with
3673
regard to accuracy and feasibility. Therefore, the OK method can be recommended for the regionalisation of substance deposition obtained by existing monitoring networks. In comparison to IDW, OK provides more accurate interpolations. This results from the crossvalidation analysis. In contrast to EDK, the OK method is less complex and needs no additional information. As a result, mapping of regional-scale distributions of the bulk deposition load of sulphate, nitrate, ammonium, sodium and the heavy metals, lead and cadmium, could be processed on the basis of the OK method. The advancement of CTMs and the opportunity to include new sources of information may allow higher resolutions in the simulation of deposition patterns and process modelling in the future. But currently, the investigated CTMs cannot be recommended in view of the feasibility for water resources management issues.
Acknowledgements We thank the La¨nderarbeitsgemeinschaft Wasser (LAWA) for funding this study. We are also grateful for uncomplicated data supply by the Saxony State Office for Forests, the Saxony State Office for Environment and Geology, the Institute of Plant and Wood Chemistry of the Dresden University of Technology, the Institute for Tropospherical Research e.V. Leipzig, the German Environmental Agency, the State Office for Ecology of Lower Saxony and the Forest Research Office of Lower Saxony. We are grateful for numerous discussions and comments of experts. Special thanks to Thomas Gauger (University of Stuttgart) and Uwe Haberlandt (formerly of the Potsdam Institute for Climate Impact Research).
References Ahmed, S., De Marsily, G., 1987. Comparison of geostatistical methods for estimating transmissivity using data on transmissivity and specific capacity. Water Resources Research 23 (9), 1717–1737. Bleeker, A., Draaijers, G.P.J., Klap, J.M., van Jaarsveld, J.A., 2000. Deposition of acidifying components and base cations in Germany in the period 1987–1995. National Institute of Public Health and Environment Protection (RIVM), Bilthoven Report no. 722108027. Davis, B.M., 1987. Uses and abuses of cross-validation in geostatistics. Math. Geological 17, 563–586. Deutsch, C.V., Journel, A.G., 1998. GSLIB—Geostatistical Software Library and User’s Guide, second ed. Oxford University Press, Oxford. Draaijers, G.P., van Leeuwen, E.P., de Jong, P.G.H., Erisman, J.W., 1996. Deposition of base-cations in Europe and its role in acid neutralisation and forest nutrition. National
ARTICLE IN PRESS 3674
F. Reinstorf et al. / Atmospheric Environment 39 (2005) 3661–3674
Institute of Public Health and Environment Protection (RIVM), Bilthoven Report no. 722108017. Draaijers, G.P., Erisman, J.W., Spranger, T., Vel, E., 1998. Quality and uncertainty aspects of forest deposition estimation using throughfall, stemflow and precipitation measurements. TNO-report TNO-MEP—R 98/093, Apeldoorn, 73pp. Englund, E., Sparks, A., 1991. GEO-EAS 1.2.1—Geostatistical Environmental Assessment Software. User’s Guide, Environmental Monitoring Systems Laboratory Office of Research and Development US Environmental Protection Agency Las Vegas, Nevada 89119. Erisman, J.W., Draaijers, G.P.J., 1995. Atmospheric deposition in relation to acidification and eutrophication. Studies in Environmental Science 63. Elsevier Science B.V., Amsterdam. Gauger, Th., Ko¨ble, R., Smiatek, G., 1997. Kartierung kritischer Belastungskonzentrationen und—raten fu¨r empfindliche O¨kosysteme in der BRD und anderen ECELa¨ndern. Endbericht zum Forschungsvorhaben 106 01 061 Teil 1: Deposition Loads. Institut fu¨r Navigation der Universita¨t Stuttgart, im Auftrag des UBA. Gauger, T.h., Ko¨ble, R., Anshelm, F., 2000. Kritische Luftschadstoff-Konzentrationen und Eintragsraten sowie ihre U¨berschreitung fu¨r Wald und Agraro¨kosysteme sowie naturnahe waldfreie O¨kosysteme. Teil 1: Deposition Loads 1987–1989 und 1993–1995. Institut fu¨r Navigation der Universita¨t Stuttgart, Endbericht 297 85 079 im Auftrag des UBA. Gauger, Th., Anshelm, F., Schuster, H., Erisman, J., Vermeulen, A., Draaijers, G., Bleeker, A., Nagel, H.-D., 2002. Mapping of ecosystem specific long-term trends in deposition loads and concentrations of air pollutants in Germany and their comparison with Critical Loads and Critical Levels. Final Report 299 42 210. Institute for Navigation, University Stuttgart. Haberlandt, U., Kite, G.W., 1998. Estimation of daily spacetime precipitation series for macro scale hydrological modelling. Hydrology Process 12, 1419–1432. Ha¨rdle, W., Simar, L., 2003. Applied Multivariate Statistical Analysis. Springer, Berlin Heidelberg New York. Hudson, G., Wackernagel, H., 1994. Mapping temperature using kriging with external drift: theory and an example from Scotland. International Journal of Climatology 14, 77–91. Isaaks, E.H., Srivastava, R.M., 1989. An Introduction to Applied Geostatistics. Oxford University Press, Oxford.
Ko¨ble, R., Smiatek, G., Gauger, Th., 1997. Kartierung kritischer Belastungskonzentrationen und—raten fu¨r empfindliche O¨kosysteme in der BRD und anderen ECELa¨ndern. Endbericht zum Forschungsvorhaben 106 01 061 Teil 2: Critical Levels. Institut fu¨r Navigation der Universita¨t Stuttgart (1997), im Auftrag des UBA. Knoth, O., Wolke, R., 1998. An explicit-implicit numerical approach for atmospheric chemistry-transport-modelling. Atmospheric Environment 32, 1785–1797. Lehmann, W., 1995. Anwendung geostatistischer Verfahren auf die Bodenfeuchte in la¨ndlichen Einzugsgebieten. IHW Institut fu¨r Hydrologie und Wasserwirtschaft Universita¨t Karlsruhe Heft 52. Matheron, G., 1971. The Theory of Regionalized Variables and its Applications. Les Cahiers du Centre de Morphologie Mathematique, Fasc.5. OMKAS, 2000. Air Hygienic Situation and Forest Condition in the Black Triangle. Proceedings of the Final Project Meeting, Saxon State Authority for Environment and Geology. Schlu¨nzen, K.H., Bigalke, K., Lenz, C.-J., Lu¨pkes, C., Niemeyer, U., von Salzen, K., 1996. Concept and Realization of the Mesoscale Transport and Fluid Model METRAS. METRAS, Technical Report 5. Meteorological Institute of the University Hamburg. Sinowski, W., 1995. Die dreidimensionale Variabilita¨t von Bodeneigenschaften—AusmaX, Ursache und Interpolation. Dissertation an der TU Mu¨nchen. Stockwell, W.R., Kley, D., 1994. The Euro-RADM mechanism: A gas-phase chemical mechanism for European air quality studies. Berichte des Forschungszentrums Ju¨lich, 2868. UBA, 1995. Monitoring der Schwermetallbelastung in der Bundesrepublik Deutschland mit Hilfe von Moosanalysen. Texte 31/95 UBA. UN ECE, 1996. Manual on methodologies for Mapping Critical Loads/Levels and geographical areas where they are exceeded. Umweltbundesamt Texte 71/96, Berlin. Van Jaarsveld, J.A., 1995. Modelling of the long-term atmospheric behaviour of pollutants on various spatial scales. National Institute of Public Health and Environment Protection (RIVM), Bilthoven Report no. 722501005. Van Leeuwen, E.P., Draaijers, G.P.J., Erisman, J.W., 1996. Mapping wet deposition of acidifying components and base cations over Europe using measurements. Atmospheric Environment 30 (14), 2495–2511.