Changes in continental and sea-salt atmospheric loadings in central Greenland during the most recent deglaciation: model-based estimates

Abstract

By fitting a very simple atmospheric impurity model to high-resolution data on ice accumulation and contaminant fluxes in the GISP2 ice core, we have estimated changes in the atmospheric concentrations of soluble major ions, insoluble particulates and ¹⁰Be during the transition from glacial to Holocene conditions. For many species, changes in concentration in the ice typically overestimate atmospheric changes, and changes in flux to the ice typically underestimate atmospheric changes, because times of increased atmospheric contaminant loading are also times of reduced snowfall. The model interpolates between the flux and concentration records by explicitly allowing for wet- and dry-deposition processes. Compared to the warm Preboreal that followed, we estimate that the atmosphere over Greenland sampled by snow accumulated during the Younger Dryas cold event contained on average four–seven times the insoluble particulates and nearly seven times the soluble calcium derived from continental sources, and about three times the sea salt but only slightly more cosmogenic ¹⁰Be.

Introduction

One goal of ice-core studies is to estimate past atmospheric concentrations of certain soluble and insoluble contaminant species. Concentrations of these species can be measured in ice cores and fluxes to the ice sheet can be calculated in some cases. However, the relation between and either concentration in the ice or flux to the ice is not direct (Davidson, 1989).

Some contaminant flux to an ice sheet occurs by wet deposition (the contaminant falls within or attached to a snowflake, either because the contaminant served as a condensation nucleus or because the contaminant was scavenged by the falling snowflake) and some by dry deposition (air-to-snow transfer without an associated water transfer; e.g. Davidson, 1989). As discussed in the Model section below, dry deposition of a contaminant species increases with its atmospheric concentration; all other things being equal, more material will fall out of air more impurities. If only dry deposition occurred, the flux of some contaminant to the snow surface would be proportional to its atmospheric concentration.

Wet deposition increases with atmospheric concentration and with the ice-accumulation rate; air with more impurities will produce snowflakes with more and more snowflakes will bring down more impurities. If only wet deposition occurred, the concentration of a contaminant in snow would be proportional to its concentration in the atmosphere sampled by that snow. In the general case, with both wet and dry deposition, the balance between them must be assessed to derive the atmospheric concentration from measurements on the ice.

Here, we estimate variations in atmospheric concentrations over millennia by extending a simple model for wet and dry deposition. We test this model using ¹⁰Be data and then apply the model to estimate changes in atmospheric loadings of sea salt and continentally derived materials through the Younger Dryas oscillation.

Data

The Greenland Ice Sheet Project II (GISP2) deep ice core was collected 28 km west of the summit of the Greenland ice sheet (72.6° N, 38.5° W, 3200 m elevation (Hodge and others, 1990)). Modern mean annual temperature is about −31°C, modern mean accumulation is about 0.24 m year⁻¹ of ice and melting occurs only about once per century or less frequently (Alley and others, 1993; Alley and Anandakrishnan, 1995).

Many of the data used here have been presented and discussed elsewhere (Alley and others, 1993; Mayewski and others, 1993a) The ¹⁰Be data are presented for the first time and will be expanded upon in subsequent publications (manuscripts in preparation by R. C. Finkel, K. Nishiizumi, and others). The particulate data were summarized briefly by Mershon and Zielinski (1993) and will be published elsewhere (M.Sc. thesis in preparation by G R. Mershon; paper in preparation by G.R. Mershon and G.A. Zielinski). We require contaminant concentrations in the ice (in units such as molecules or particles m ⁻³), fluxes to the ice (molecules or particles m⁻² year⁻¹, calculated as the product of the concentration and the ice accumulation or flux in m year⁻¹) and the age relative to the major late-glacial climatic events (e.g. Johnsen and others, 1992; Alley and others, 1993; Taylor and others, 1993b).

Ice-core ages were determined by counting annual indicators (visible strata, oscillations in electrical conductivity, variations in laser-light scattering from dust, plus oxygen-isotopic and chemical variations at certain depths). Absolute-age estimates were checked against fall-out from historically dated volcanoes over the most recent 2000 years and against independent ages for major events during the deglaciation. Assigned ages are believed to be accurate to better than 3% through the Late Glacial and 1–2% through the Holacene (Taylor and others, 1992; Alley and others, 1993).

Accumulation rates were estimated Gom annual-layer thicknesses by correcting for the layer thinning caused by the spreading and stretching of the ice sheet during its flow. Layer thicknesses predicted for the ice core by Schøtt and others (1992) were exceedingly accurate, so we simply scaled the flow corrections from their model (Alley and others, 1993). Further work on ice-flow corrections (Bolzan and others, 1995; Cutler and others, 1995) in general supports the Schøtt and others (1992) reconstructions. The time intervals considered here are short enough that any “drift” between different reconstructions owing to ice-sheet thickness changes should be small enough to be ignored. All data discussed here are from much shallower than the zone of flow disturbances near the bed (Taylor and others, 1993a; Alley and others, 1995).

Data on concentrations of soluble major anions and cations (sulfate, nitrate, ammonium, chloride, sodium, potassium, magnesium and calcium) are from Mayewski and others (1993a), who interpreted the high-resolution time series of these species. Approximately bi-yearly samples were analyzed, ranging from 0.07 to 0.21 m in length and containing one to a few years but typically 2–3 years.

Insoluble particulates were measured in a class 100 clean room using an Elzone 280PC (Coulter-type) particle counter equipped with a 30.308 μm orifice. Grain-size determinations were made in 64 separate channels ranging from 0.67 to 12.88 μm. The 64 channels can be grouped into five bands without much loss of detail, avoiding some noise associated with the narrow channels and also simplifying data analysis and presentation (Mershon and Zielinski, 1993; paper in preparation by G. R. Mershon). We briefly summarize from those five bands and from the finest channel within the finest band (0.67–0.70 μm); full results will be presented elsewhere (paper in preparation by G.R. Mershon and G, A. Zielinski). Major-ion and particulate sampling used the same depth ranges and sample sizes. The eluent used for the chemical samples was a weak acid that would have increased dissolution of calcium carbonate dust; thus, some dust should be measured in both data sets.

¹⁰Be was measured on samples of core melted in the field. After the addition of carrier, Be was concentrated from 1 to 2 kg of meltwater using ion-exchange chromatography. The Be was purified at the University of California–San Diego and the ¹⁰Be measured by accelerator mass spectrometry at the Lawrence Livermore National Laboratory Center for Accelerator Mass Spectrometry (Davis and others, 1990). Samples ranged from 1.4 to 3.4 m in length or typically 25 to 100 years per sample.

Major-ion and particulate data extend from 11 322 to 14035years before AD 1950 (ABP), or 1657–1762 m depth in the ice core, which spans the interesting time period of the latter part of the Bolling/Alleröd (BA) warm event, the entire Younger Dryas (YD) cold event and the earlier part of the warm Preboreal (PB) that followed the Younger Dryas (Fig. 1).

Fig. 1. Time series of fluxes of calcium (top panel), chloride (second down), ¹⁰Be (third down) and water (ice accumulation; bottom panel) plotted against age in years before AD 1950 for the GISP2 ice core. Concentrations in each data set have been normalized by dividing by the mean of that data set to allow display on similar scales. We extended the ¹⁰Be and ice-flux records somewhat further than the chloride and calcium records to test the effect of the extra data points on the statistical validity of the results, as described in the text. The accumulation-rate record has-been smoothed to the sampling internal of the chemical species to which if was compared; for the younger samples, we display the resolution of the calcium and chloride samples; for the older part, we show the lower resolution of the ¹⁰Be samples.

To consider separate climate states, we sub-set the particulate and chemical data and eliminated regions near transitions (11 637–11 759 ABP or 1675–1680 m for the PB/YD transition and 12 854–13 049 ABP or 1710–1718 m for the YD/ΒΑ transition) to leave samples clearly representative of the time periods under consideration. This leaves sub-sets of 145 PB samples, 299 YD samples and 232 BA samples for which we have major-ion chemical data, particulate data and average ice fluxes. For ¹⁰Be, the major effect of eliminating data near the transitions is to reduce the statistical confidence without changing the results, so we retained all of the data. We also examine some ¹⁰Be data extending through the BA and into the Glacial Maximum or Oldest Dryas (OD) (Fig. 1).

Model

Model principles

The simplest plausible model (e.g. Legrand, 1987; Davidson, 1989; Whung and others, 1994) relating contaminant flux to atmospheric concentration and ice flux is

(1)

Here, f _i, is the flux of contaminant species i, C _i, is the atmospheric concentration of this species, b is the ice-accumulation rate, k ₁ is the dry-deposition velocity, k ₂ is the dimensionless scavenging ratio for precipitation and K_i is the concentration of the contaminant species in the ice core.

Assume initially that k ₁ and k ₂ are constants (we discuss this assumption in the Model test section, below). Then Equation (1) shows that, if dry deposition dominates (k ₂ C_ib ≪ k ₁ C _i), the flux of some contaminant to the ice sheet will be proportional to its atmospheric concentration. If wet deposition dominates (k ₂ C_ib » k ₁ C _i), the concentration of some contaminant in the ice sheet will be proportional to its concentration in the atmosphere sampled by the accumulated snow. In the general case with both wet and dry deposition, an increase in snow accumulation with constant atmospheric concentration of a contaminant will increase the flux of that contaminant to the ice sheet but decrease the concentration of that contaminant in the accumulated snow.

Data on ice-core contaminants usually are reported as concentrations or, where snow accumulation is known accurately enough, as concentrations and fluxes (e.g. Mayewski and others, 1993a, b). For sites dominated by wet deposition, the time series of concentration provides the better estimate of atmospheric loading, whereas for sites dominated by dry deposition the time series of flux provides the better estimate. Here, we present the simplest physically based method to interpolate between changes in flux and changes in concentration so as to estimate changes in contaminant atmospheric loading for sites with both wet and dry deposition.

During a climatic regime such as the Younger Dryas, the atmospheric concentration will have some average value, C _i . Ice accumulation will vary from year to year, owing in part to local effects. Years of locally high ice accumulation will have a larger contaminant flux than years of locally low ice accumulation, because wet deposition will bring down more contaminants during the high-accumulation years than during the low-accumulation years. Thus, if k ₁ and k ₂ are constants (see below), a plot of contaminant flux (ƒ _i) vs ice accumulation (b) should produce a straight line with slope k ₂ C _i and intercept k ₁ C _i The intercept (the contaminant flux extrapolated to zero ice accumulation) can be taken as an approximation of the dry-deposition rate (see Discussion) and the additional flux is wet deposition.

We expect deviations of individual data points from this line owing to variations in the atmospheric concentration, C_i , from its mean value, C _i . Deviations may also be caused by measurement error, time variation of k ₁ or k ₂, or other processes not included in the model. Errors in chemical-concentration and particle-concentration measurements are well-characterized and not large. Misidentification of annual markers in the ice core is possible. However, in recent times when we can check our main counting methods against independent annual-layer indicators and absolute-time horizons, our counts are quite accurate (Taylor and others, 1992, Alley and others, 1993); thus, we believe that errors from this source are not large or at least are not systematic.

Model applicability

The model in Equation (1) clearly is not as complex as However, this model provides a good fit to data from southern Greenland on methanesulfonate fluxes over the last few centuries (Whung and others, 1994). It also accurately fits modern data on spatial distribution of fluxes of sea salt and sulfate, and to a lesser extent nitrate, in East Antarctica (Legrand, 1987).

Ice-flux/contaminant-flux (b, f) plots of our data for different climatic events and for different contaminants generally show positive slopes, as expected from the model. For each of eight soluble ionic species and six insoluble particulate size classes in each of three time periods (the Bolling/Alleröd (BA), Younger Dryas (YD) and Preboreal (PB) sub-sets as described above), simple regression analysis on (b, f) plots produces positive slopes for all of the PB and YD cases and 12 of the 14 BA cases, with all but one of the PB and YD cases different from zero with >95% confidence based on the standard t-test. (Only five of the BA cases are significantly different from see below.) Because f is calculated as the product of concentration, K, and accumulation rate, b, a tendency for positive slope is built into the analysis. However, we can exclude many possible circumstances that would produce zero or negative slopes and invalidate our model.

For example, if wet deposition were extremely inefficient, contaminant flux would be independent of ice flux, which is not observed. Similarly, if wet deposition were extremely efficient, then contaminant flux would become independent of ice flux with increasing ice flux, yielding a horizontal line for high b on a (b, f) plot. We observe no such tendency for curvature to horizontal with increasing b. Mixing of samples from different populations could yield (b, f) data arrays with negative slopes if contaminant source strengths are higher for populations with lower average ice flux. When we combine the PB, YD and BA populations into a single group, we do obtain negative slopes for most contaminants. For separate times, the BA accounts for all except one of the (B, f) slopes that are not positive with >95% confidence, suggesting mixing of populations in the BA only. Indeed, inspection of the data sets shows that much of the BA is similar to the PB but with short YD-like events, also suggesting mixing of populations.

Correlation coefficients, r, for the (b, f) linear regressions for different times and contaminants range from 0.6 down to 0.1 and average 0.3, thus, significant variability not included in our model must be present, as expected. The reader will recall from the standard t-test (e.g. Till, 1974), the minimum significant correlation coefficient for data sets with a large number of points, n, scales approximately with 1/√n. For our data sets and for a one-tailed test with 95% confidence, r >≈ 0.1 is significant.

To compare climatic regimes from different time intervals t = u and t = v, we assume that the processes of contaminant removal from the atmosphere did not change and thus that k ₁ and k ₂ are constants. (We test this assumption in the next section and find that it is consistent with available data for one species.) Then the ratio, Ri^u:v , of the slopes on a (b, f) plot from the two regimes should equal the ratio of the intercepts and of the average atmospheric concentrations

(2)

The two lines for the two time intervals have three rather than four free parameters: one slope, one intercept and the ratio Ri^u:v . Of these three parameters, the ratio is most directly interpretable and so of greatest paleoclimatic interest.

To conduct the joint linear regression with the constraint in Equation (2), we use a simple inverse technique based on a generalization of Newton’s method to minimize the total variance of the model. The diagonal elements of the covariance matrix returned by the inversion are the variances of the model parameters (see Press and others, 1988). We find that the lines produced this joint linear regression are statistically indistinguishable from those produced by separate linear regressions from the two time periods considered; examples are included in Figure 2 and are listed in Tables 1 and 2. We favor the joint regression based on physical reasoning and on the reduction in uncertainty gained from this physical reasoning but we tabulate some results of separate regressions for readers who may be more familiar with this treatment. As with the individual regressions, we find statistically significant results using this joint-regression technique.

Fig. 2. Normalized fluxes of ¹⁰Be (a), calcium (b), and chloride (c), plotted against normalized flux of water (ice-accumulation rate) for approximately bi-yearly samples. Populations representing the Preboreal (PB) and Younger Dryas (YD) are shown for calcium, and chloride. Bolling/Alleröd (ΒA) and Oldest Dryas (OD) points also, art shown for ¹⁰Be; ΒA data plot between OD and YD for calcium and chloride. Lines obtained by the joint regression of YD and PB data using our model are shown dashed extending to their intercepts for PB and YD for calcium and chloride, and for PB + BA and YD + OD for ¹⁰Be. The shorter, solid-line segments in (b) and (c) are the individual regression lines for the YD only and for the PB only.

Table 1. Ratios between Younger Dryas (YD), Preboreal (PB) and Bolling/Alleröd (BA) atmospheric chemical concentrations, together with the standard deviations of the ratios, determined from the joint regression with the ratio of the slope equal to the ratio the intercepts. Because ratios were calculated between pairs of times (YD: PB, YD: ΒA, ΒA: PB) rather than fitting all three times at once, the product (YD : BA) × (Β A : PB) is not identical to the single ratio (YD : PB), although they are similar. We also calculated, for ¹⁰Be only, YD: (PB + BA) = 1.21 ± 0.89 and (YD + OD): (PB + BA) = 1.27 ± 0.72 where OD is the cold Oldest Dryas; in Figure 3 we use physical as well as statistical constraints to obtain (YD + OD) : (PB + BA) = l.22±0.39. The main effect of using these expanded cold: warm data sets for ¹⁰Be is to reduce the statistical uncertainty of the ratios. The ratios of slopes from separate regressions are also shown for the YD: PB. case; for ¹⁰Be YD + OD : PB + BA. this yields 1.52±0.70. The separate regressions tend to yield slightly larger ratios than the joint regression but in no case is the difference highly significant

Fig. 3. Alternative way to estimate R^{coulD: warm} _10Be. The mean values of normalized chemical and ice fluxes for warm (PB + ΒA) and for cold (YD + OD) periods, and their uncertainties, are shown by the small crosses. Regression lines through these means, having the ratio of their slopes equal to the ratio of their intercepts and having positive slopes and intercepts, are limited by the solid and dashed lines shown. This range, when increased to allow for the small uncertainty in the means, yields R^{coulD: warm} _10Be = 1.22 ± 0.39.

Table 2. Ratios between Younger Dryas (YD), Preboreal (PB) and Boiling/Άlleröd (BA) atmospheric concentrations of insoluble particulates, together with the standard deviations of the ratios, determined from the joint regression with the ratio of the slopes equal to the ratio of the intercepts. Because ratios were calculated between pairs of times (YD : PB, YD: ΒA, BA: PB) rather than fitting all three times at once, the product (YD : ΒΑ) × (BA : PB) is not identical to the single ratio (YD: PB), although they me similar. The ratios of slopes from separate regressions are also shown for the YD : PB case. The separate regressions tend to yield slightly larger ratios than the joint regression but in no case is the difference highly significant

Model test

We first test the model against ¹⁰Be data and show that any variations in k ₁ and k ₂ between the cold YD and the warm PB are sufficiently “small” to be ignored for this species. In subsequent sections, we use this result and assumptions about behavior of other species to estimate that C changes about three-fold for sea salt, seven-fold for continentally derived soluble calcium and four-fold to seven-fold for insoluble particulates from YD to PB; thus, “small” variations in the constants are those less than a few tens of per cent. We use ¹⁰Βe for the model test, because we can estimate the ratio R _10Be ^{YD: PB} using independent means, to test the model parameters derived from GISP2 ¹⁰Be data.

¹⁰Be is a long-lived (1.5 × 10⁶year) radioactive nuclide produced by cosmic rays, mostly by spallation reactions in the atmosphere, with about two-thirds of the production in the stratosphere and one-third in the troposphere (Lal and Peters, 1967). Residence time is almost a year in the stratosphere but only about 3 weeks in the troposphere.

¹⁰Be production can change in response to changes in the galactic cosmic-ray flux (for example, caused by supernovas), or in response to changes in shielding of the Earth from cosmic rays (caused by changes in solar activity or in the Earth’s magnetic field). Small variations in production occur correlated with sun-spot activity (Beer and others, 1983, 1990), events of increased production have occurred about 35 000 and perhaps 60 000 ABP, probably in response to supernovas or magnetic field events (Raisbeck and others, 1987, 1992) and slow variations in production similar to those for ¹⁴C may have occurred in response to slow changes in the Earth’s magnetic field (Lal, 1992). However, production probably has been nearly constant for our purposes (that is, variations on the order of 10% or less) over the most recent deglaciation including our study interval (McHargue and Damon, 1991; Lal, 1992; Mazaud and others, 1992; Raisbeck and others, 1992).

Figure 1 shows the time series of ¹⁰Be flux during the deglaciation and Figure 2a shows how ¹⁰Br varies with accumulation rate for the various time Intervals considered. It appears that both dry deposition (the intercept in Figure 2a) and wet deposition are important in the ¹⁰Be flux in central Greenland.

Because of the likelihood that ¹⁰Be production has been nearly constant, we can use Figure 2 to test whether our model fits the data, and thus whether k ₁ and k ₂ have been (nearly) constant over time. We do so by using simple arguments about atmospheric processes to show that changes in the behavior of the ¹⁰Be flux between cold (glacial or stadial) and warm (interglacial or interstadial) times must fall between rather narrow limits if k ₁ and k ₂ are constants and that the observed behavior falls between these limits.

Over most of the world (where precipitation rates are higher than in central Greenland), most ¹⁰Be is removed from the atmosphere by wet deposition (McHargue and Damon, 1991). During cold periods, such as the Younger Dryas or the Oldest Dryas, results from atmospheric general circulation models indicate that the global precipitation rate was either unchanged or reduced by only a few per cent from warm periods such as today (e.g. Kutzbach and Guetter, 1986; Lautenschlager and Herterich, 1990; Joussaume, 1993). If ¹⁰ Be were globally mixed in the troposphere, its atmospheric concentration would be similar during cold and warm periods because production nor removal rates would have changed much over time. During cold periods, the locally reduced snowfall in central Greenland would cause the ¹⁰Be flux to be reduced there. (Because Greenland receives such a small fraction of the global flux, a change in deposition in Greenland would have a negligible effect on the atmospheric loading of a globally mixed species.) Assuming no change of k ₁ and k ₂, we would End in this constant-atmospheric-concentration limit that R ^{coulD: warm} _10Be (const, cone.) = 1.

However, the 3 week tropospheric residence time of ¹⁰Be is not long enough for complete global mixing in the troposphere. Suppose instead that no lateral mixing of ¹⁰Be occurs in the atmosphere, so that the flux of ¹⁰Be to Greenland is constant over time. If wet deposition is significant, the reduced precipitation and wet deposition during cold periods would increase the atmospheric residence time and thus the atmospheric concentration of ¹⁰Be in air sampled by Greenlandic accumulation. The increased atmospheric concentration would in crease the dry-deposition flux. Observing from Figure 2a that dry and wet deposition are about equal in cold periods, we can calculate from Equations (1) and (2)th at this constant-flux limit with k ¹ and k ₂ constant would yield R ^{coulD: warm} _10Be (const. flux.) = l.5 between cold and warm periods. (If all flux were dry deposition at all times, changes in precipitation from cold to warm periods would not affect the atmospheric concentration, giving R = 1 between cold and warm times. If all flux were wet deposition at all times, then constant flux with halved precipitation would require doubled concentration in cold times, or R = 2, assuming a linear model. With both dry and wet deposition important, R flls between these. We make the assumption that precipitation and accumulation are approximately equal here.)

Again, this limit is unlikely to be correct. The modern tropospheric residence time of about 3 weeks (Lal and Peters, 1967), which may have lengthened during the Younger Dryas, is enough to allow mixing to occur beyond the immediate North Atlantic Ocean, where the Younger Dryas cold event is believed to have been strongest (e.g. Rind and others, 1986). Thus, we expect that 1 < R ^{coulD: warm} _10Be < 1.5 for ¹⁰Be between the Younger Dryas and Preboreal times, if k ₁ and k ₂ are constants.

Thus far, our data set includes many ¹⁰Be samples from the cold Younger Dryas (n = 21), the warm Boiling/Alleröd (n = 28) and the cold Oldest Dryas (n = 22) but only a few from the warm Preboreal (n = 5). We have applied our model to various combinations of these, as shown in Table 1. The result is that 1 < R ^{coulD: warm} _10Be < 1.5, as expected, with a suggestion that R falls towards the lower end of this range.

To improve the statistical confidence of this result, we also have conducted the analysis while enforcing the physical constraint that the dry deposition and wet deposition are positive (i.e. slopes and intercepts are positive; Fig. 3). We calculate the mean and the uncertainty (s/√n, where s is the sample standard deviation; e.g. Till, 1974, chapter 4) of cold-period and warm-period data on the chemical-flux/ice-flux plot. If we construct pairs of lines with positive slopes and positive intercepts through these means such that the ratio of their slopes equals the ratio of their intercepts, the limiting values for R ^{coulD: warm} _10Be occur as dry deposition (the intercept) approaches zero and as wet deposition (the slope) approaches zero. If we then calculate R for these limiting cases, and allow For the (s/√n uncertainty in the mean values, we obtain R ^{coulD: warm} _10Be = 1.22 ± 0.39. This is almost identical to our previous result but with smaller uncertainty because of greater constraint by physical reasoning.

The observed behavior of ¹⁰Be and ice flux from cold to warm periods (R = 1.22 ± 0.39) is statistically indistinguishable from the behavior expected if our model is correct (R = 1.25 ± 0.25). We conclude that any variations that have occurred in k ₁ and k ₂ for ¹⁰Be have been small, on the order of tens of per cent or less.

The tendency for R ^{coulD: warm} _10Be to be at the low end of its possible range is consistent with air masses affecting Greenland that experienced smaller reductions in precipitation during cold periods than did Greenland; however, this result is statistically valid at less than the one standard-deviation level, so we cannot insist on it. Further, more-careful analysis of the ¹⁰Be data is in progress and may allow more-refined conclusions to be drawn.

Results

Insoluble particulates

We now wish to apply this model to other species. However, some care is required. The ¹⁰Be used in our model test is believed to be deposited attached to “accumulation-mode” (0.1–1 µm) sulfate particles, which probably have had nearly constant grain-size distribution through the entire period of interest (reviewed by McHargue and Damon (1991)). Strictly, we have tested the model only for these particles. Through a string of explicit assumptions (which we argue are reasonable assumptions), we now apply our model to other sizes and types of particles and to bulk chemistry. The analysis path is sketched in Table 3.

Table 3. Path followed to achieve conclusions of this paper. The model tests and assumptions outlined here are described more fully in the text

Results of application of the model to the data are shown in Table 2 for five size classes of insoluble particulates and for the smallest particulate channel measured (0.67–0.7 µm). The smallest channel, and the smallest of the five size classes, fall within the accumulation mode presumably occupied by the ¹⁰Be-bearing sulfate particles. If we assume that the model is valid for particles in the accumulation mode, we find that the atmosphere sampled by snow acrumulated during the Younger Dryas contained about four times the insoluble particulates of the Preboreal in this size range. If we assume that the model is valid for comparisons within any size class, we find that the change in atmospheric loading from Younger Dryas to Preboreal was larger for coarser particles than for finer particles and exceeded seven-fold for the coarsest particles measured here.

The data in Table 2 show that the insoluble particulates in the atmosphere were on average coarser during the YD than during the PB. This is consistent with a variety of data comparing particulates in ice from glacial and interglacial climates in Greenland and Antarctic ice cores (e.g. Thompson, 1977; Petit and others, 1981; Menhon and Zielinski, 1993). The model parameters k 1 and k 2 probably are increasing functions of grain-size for particle sizes that are significant in mass loading of the atmosphere (e.g. Warneck, 1988; Schumann, 1991; Hillamo and others, 1993). Thus, flux of contaminants to the ice, and concentration in the ice, may increase from to cold times because of increased atmospheric removal efficiency of coarsened grain-sizes during cold times as well as because of higher atmospheric loadings during cold times.

Because insoluble particulates are measured in hands with non-zero width, coarsening of the grain-size distribution shifts the average size within an interval slightly. However, we calculate that this effect is not significant in our results, except possibly for the coarsest interval, which is broader than the others. Most measures of grain-size distribution for atmospheric contaminants form smooth curves when plotted against the logarithm of particle radius. We plotted the results for the five broad bands shown in Table 2 against the particle radius, drew a smooth curve through them and then compared this smooth curve to results from analysis of narrower single channels.

The general result is that the broad bands do slightly overesimate the YD: PB change, as expected. For the finest channel, this overestimate amounts to about 2% of the ratio (3.96 vs ≈4.02). For various channels in the coarsest band, the overestimate averages about 10%. In comparison, the curve-fitting errors are typically 20% of the ratios, so these differences are not highly significant. The coarsest band is the widest considered and so would be expected to have the largest errors. The coarser individual channels have very few particles in them per sample, causing the curve-fitting errors to be larger than for the bands combining several channels.

Next, we wish to apply our model to data on soluble contaminants. However, our measurements cannot yield particle-size distributions for soluble species but only total concentrations in the ice. This is equivalent to having only a single band for particulates. (We are confident that the particulate-size classes sampled capture the particulates that contribute significantly to the mass flux.) Again, if the gram-size distribution of the particles carrying soluble species coarsened between climate states and the coarser particles are removed more efficiently, our model would misinterpret the increased removal of the coarser particles as an increased atmospheric loading. The atmospheric sampling efforts under way in central Greenland (eg. Bergin and others, 1994; Dibb and others, 1994) eventually should allow model-based assessments of the magnitude of this effect. Pending further results from such studies, we cannot fully account for the effects of grain-size changes but we can use our particulate data to estimate the magnitude of possible errors.

To treat our insoluble-particulate data as if they were collected in the same manner as soluble-chemical data, for each sample we summed the particulates in each channel to obtain the total mass and then we applied the model to these summed data. Comparing the YD to the PB, we obtain a ratio of 6.47 for these particulate-mass data.

If we knew the size distribution of particles in the atmosphere at any time, we could then calculate the change in atmospheric mass loading and the bias introduced by gram-size-dependent removal processes. For example, we find that the 0.67–0.70 μm particulate channel shows a YD : PB ratio of 3.96. If the atmospheric mass loading during YD and PB were almost entirely of particles in this size class (with coarser particles being important in the ice only because of an extremely strong size-dependence of k ₁ and k ₂), then our bulk estimate of a YD : PB change in atmospheric mass loading of 6.47 would be about 60% higher than the actual change of 3.96. If 1 μm particles dominated the atmospheric mass loading, the calculated 6.47 would exceed the actual 4.33 by almost 50%.

However, most data suggest that the atmospheric loading of continentally derived material and sea salt is dominated by slightly coarser particles. Hillamo and others (1993) showed that for sea-salt and continental aerosols at Dye 3, southern Greenland, sampled during March 1989, the mass distribution typically showed a peak near 2 μm. This is consistent with a range of other results from remote sites (see Warneck, 1988, chapter 7). If we note that grain-size distributions in the ice core from the PB are slightly coarser than those from recent times, and those from the YD are even coarser, it is reasonable to suggest that the change we calculate in atmospheric loading for 2 µm particles or slightly coarser particles comes close to the actual change in mass loading in the atmosphere. The difference between the change m this size and the change calculated by application of our model to total particle masses in the ice is the result of particle-size-dependent removal processes combined with changes in particle-size distribution.

Two channels come closest to 2 μm size: the 1.4–2.0 μm channel, with a YD : PB ratio of 5.24 and the 2.0–3.0 µm channel, with a YD: PB ratio of 5.86, vs the model-obtained result for total particulate massses in the ice of 6.47 By fitting smoothed curves to number-concentration data plotted against the logarithm of the particle sizes considered (usually referred to as dN/dlnR; e.g. Warneck, 1988), we obtain a YD : PB ratio of 5.5 for the atmospheric loading of 2 μm particles. This is about 17% lower than we obtain from bulk analysis of the particles. We therefore suggest that analysis of the bulk data overestimates the change from PB to YD by approximately this much or less.

Several factors may play a role in the small magnitude of this effect. If processes such as filtration of particle-laden air below the snow surface (Cunningham and Waddington, 1993; Hillamo and others, 1993) or scavenging by falling snowflakes are sufficiently vigorous, they may remove essentially all particles from the air so that size-dependent fractionation is unimportant. The change in grain-size of particles appears small (only a few tenths of a micron between the mass-weighted mean sue in the ice for YD and PB) and may not be large enough to affect the results significantly.

The change in particle-size distributions from cold to warm periods was probably larger for the continentally derived materials that produce the insoluble particulates we study than for sea salt or most other chemicals (reviewed in Warneck (1988)). We therefore suggest that application of our model to bulk measures of contaminant loading, such as are obtained from wet chemical analyses, produces an overestimate of changes in atmospheric loading from cold to warm periods which nonetheless is fairly close to being accurate. We expect that ongoing atmospheric research will allow better estimates of this bias but that it is on the order of 10%. We do not correct for this bias.

Soluble marine and continental major ions

For clarity, we concentrate on chloride, which is dominated by marine sources, and calcium, which is dominated by continental sources (e.g. Clausen and Langway, 1989; Delmas and Legrand, 1989; Mayewski and others, 1990), although we tabulate data on other ions for interested readers (Table 1), The C1: Na weight ratios averaged over each of the time periods are close to sea-water values (2.45 ± 0.06 PB, 1.89 ± 0.03 YD, 2.21 ± 0.06 BA, ≈1.82 sea water). C1 and Na give statistirally indistinguishable results in our analyses, so we do not address the question of which is the better marine indicator (nor can we easily resolve changes in the gas-phase behavior of C1, because the marine signal is dominant). The small sea-salt contribution to sulfate has not been corrected for here. Magnesium has marine as well as continental sources and potassium may have significant biomass-burning as well as continental-dust sources (e.g Clausen and Langway, 1989; Delmas and Legrand, 1989; Mayewski and others, 1990). In general, the behavior of magnesium and potassium in our data falls somewhere between the continental calcium and the marine sodium and chloride. Nitrate, sulfate and ammonium have more-complicated sources and atmospheric chemistry (Clausen and Langway, 1989; Delmas and Legrand, 1989; Mayewski and others, 1990; Legrand and others, 1992), and we leave consideration of them for other studies; we tabulate results for these species strictly as a service to interested readers, and offer no interpretations or conclusions.

Application of our model to the data yields the regression lines shown in Figure 2, and the ratios and uncertainties listed in Table 1. The simple interpretation is that the cold YD atmosphere over Greenland sampled by accumulated snow showed a three-fold increase in the concentration of sea salt and a seven-fold increase in the concentration of continentally derived soluble calcium, compared to the warm PB atmosphere that followed. The BA was generally warm but included much variability, and atmospheric concentrations of sea salt, soluble calcium and most of the other chemicals studied here fell between the PB and YD values.

Because we do not have samples with very low accumulation rates, the intercepts in Figure 2 require long extrapolations of the regression lines and so are not statistically well-constrained (also see Discussion). Thus we are unable to draw any strong conclusions about relative importance of dry vs wet deposition for individual species. However, the best estimate of the fraction of total flux contributed by dry deposition in the YD exceeds that for the PB for every ion considered. If we average the per cent contribution of dry deposition to total flux for all species in the PB and compare to the YD, we find that dry deposition was more important in the YD than in the PB with >90% confidence (18% YD vs 11% PB dry deposition for the mean of the best estimates for the eight major ions).

Similarly, dry deposition provided a larger fraction of the total flux of insoluble particulates during cold times than during warm times (43% YD vs 26% PB dry depostion for the mean of the best estimates of the five size bands). In addition, dry deposition a somewhat more important for insoluble particulates than for soluble contaminants and dry deposition is more important for coarser particles than for finer ones (47% YD, 30% PB for the two coarser bands, vs 38% YD, 22% PB for the two finer ones).

If dry deposition were identically zero, then we would expect concentration of a contaminant in the ice core to be proportional to its atmospheric concentration The contribution of dry deposition to total flux is typically rather low, especially in the warm periods, which means that concentrations in the ice provide fairly good estimates of concentrations in the atmosphere at this site (Table 4), although the full analysis presented here should provide a better estimate Were dry deposition dominant (as is possible for some species in central regions of East Antarctica (Legrand, 1987: Davidson, 1989), then chemical flux would track atmospheric concentration more closely than would concentration in the ice.

Table 4. Ratios between Younger Dryas (TD) and Preboreal (PB) concentrations of insoluble particulates, soluble calcium and soluble chloride, determined from fluxes to the ice sheet, our model and concentrations in the ice sheet

Discussion

It is difficult to overstate the complexity of the physical and chemical processes leading to the archival of atmospheric contaminants in an ice sheet Much progress has been made in understanding these processes (e.g, Davidson and others, 1985; Dibb and others, 1992; Jaffrezo and others, 1993) but the goal of inverting an ice-core record for the full suite of atmospheric conditions that produced it remains elusive. Pending success of that effort, simplified treatments such as this one may have a place in ice-core analyses.

Analogy might be drawn to stable isotopic paleother-mometry. The list of possible factors that affect isotopic temperatures is long—changes in source region, transport path, snow drifting, sublimation, seasonality of precipitation, etc., as well as changes in temperature. Yet, in many cases, isotopic ratios provide accurate paleothermometers (e,g. Peel and others, 1988; Cuffey and others, 1992, 1994).

Some of our assumptions merit further discussion. We equate the intercept on a water-flux/impurity-flux (b–f) plot with the dry-deposition rate. However, Legrand (1987) wrote that this was “debatable”. A first obvious difference is that Legrand considered the spatial variations of impurity and water flux but ice-core interpretation is a temporal study. As shown by Peel and others (1988) for a possibly analogous case, spatial gradients of stable-isotopic ratios may be significantly different from the temporal ones of interest.

A second factor is that, because we have not identified any years with very low water flux (b, or snow accumulation), the intercept is poorly constrained because of the need for long extrapolations. We also cannot observe any tendency for curvature of the b– f data at small b. Should curvature occur owing to a greater efficiency of atmospheric snow-out at small b, then our intercept values are curve-fitting parameters that overestimate the true dry-deposition rate. (The opposite curvature, requiring greatly inefficient snow-out at low snowfall rates, is physically implausible.) Reversible deposition of species with gaseous phases (including water vapor) also would complicate the identification of the intercept as the dry-deposition rate (Bergin and others, 1994). However, the lack of any clear curvature in the data we do have, and the success of our model test, lead us to trust our interpretation of the intercept as a good estimate of the dry-deposition rate.

We wish to emphasize that this and all of our model tests are based on our data They thus are site-specific. Before techniques such as this can be used with confidence at other sites, similar tests must he conducted. We might speculate, for example, that some of our successes are related to the low importance of dry deposition at GISP2, and that results possibly could differ at sites with a higher ratio of dry to total deposition where any errors related to identification of the intercept with the dry-deposition rate would be more important.

The low correlation coefficients for regression lines on b–ƒ plots at first might alarm readers, suggesting that the model does not fit the data. In fact, the statistical analyses give great confidence that impurity flux ƒ does increase with water flux b for almost all species and particulate-size classes in the YD and the PB, and for many in the BA. Low correlation coefficients but high t values are perfectly consistent provided there is significant variability or “noise” in the system (e.g. Till, 1974). Such variability is observed in the contaminant data and is expected.

Suppose, for example, that someone had conducted a similar study in the eastern United States over the last century. It is a reasonable expectation that a large signal would have been evident from the dustbowl of the 1930s. But, because we have taken 1000 year-long climate states, we lump dustbowl-type events with non-dustbowl years, Little Ice Age type events with Medieval Warm Period type events, and so on to produce a “noisy” data set that cannot he described completely by assuming a constant source strength and atmospheric loading (see Mayewski and others, 1993b)

Given two reference points or a reference line (for example, the slope and intercept of the regression line through some selected part of the Holocene data set), our model can be applied to a single data point to produce an estimate of changes in atmospheric loading relative to the reference state. Thus, we can estimate atmospheric loadings with the same time resolution as the sampling scheme, which may be of interest in some studies. We have chosen to focus on longer intervals because we could test our model over longer times, and because we could use the many data points within a long interval to estimate uncertainties as well as averages; however, we hope to examine high-resolution results in the future.

Many other possible complicating factors may exist for this model, including correlations between accumulation-rate variations and atmospheric-concentration variations within a climatic regime. However, the model test using ¹⁰Be data indicates that model errors are small for ¹⁰Be that they affect the YD and PB populations equally (hence they largely cancel out when YD: PB ratios are formed), or that there are highly fortuitous offsetting effects. We consider the final possibility to be highly unlikely (although we obviously will be seeking further verification as more data become available). We do not worry too much about the other two, because our results are valid under both.

We cannot, of course, absolutely exclude the possibility that the model works for ¹⁰Be but not for other species. We note that we chose to test the model against ¹⁰Be, and we estimated the limiting behaviour of ¹⁰Be, before we had analyzed the ¹⁰Be data; thus, any claim for uniqueness of ¹⁰Be would require a surprising coincidence. The availability of insoluble-particulate data within size classes allows us to eliminate grain-size-dependent changes in atmospheric-deposition processes as significant complicating factors in interpreting the particulate data, and the limited available data suggest that grain-size changes have not greatly complicated applicatron of the model to soluble-impurity data. We therefore believe that the model provides a working tool for estimating changes in average atmospheric loadings of contaminants over ice sheets.

Notice that, if dry and wet deposition sample different levels of the atmosphere with different concentrations of some contaminant at a time, the similarity of the results we obtain from joint and separate regressions (Tables 1 and 2) indicates that the proportional changes in those concentrations between climate states are similar. Were the proportional changes different, the ratio of slopes returned by the separate regressions would differ from the ratio of intercepts, causing the lines from the individual regressions to differ from those of the joint regression. The errors listed in Tables 1 and 2 are large enough that we cannot exclude the possibility of some such altitude-varying changes in concentrations but we cannot demonstrate them either.

Conclusions

Studies of contaminant concentrations in and fluxes to ice cores have provided great insights into changes in sources, transport paths and effects of aerosols (e.g. Petit and others, 1981; Herron and Langway, 1985; Harvey, 1988; Thompson and others, 1988; Mayewski and others, 1993a). However, some uncertainty has always been attached to the interpretation of atmospheric loadings because of the poorly quantified role of changes in wet and dry deposition. The availability of annually resolved ice cores from regions of simple ice flow allows us to changes in snow accumulation, and thus to estimate changes in atmospheric loadings with less uncertainty than previously. The simple model of < Equations (1) and (2) is a first attempt at doing so.

We certainly have not “solved” atmospheric-removal problems. Full physical models of removal processes (e.g. Davidson, 1989) eventually should replace regression models such as ours. However, we argue that our simple model improves on the use either of contaminant concentration in ice or of contaminant flux to ice sheets to estimate changes in atmospheric loadings.

Application of our model produces estimates that, compared to the warmer Preboreal that followed, the atmosphere over central Greenland sampled by snow accumulated during the Younger Dryas cold event contained on average about three times the sea salt, seven times the continentally derived soluble calcium and four times (for sub-µm dust) to seven times (for several-μm dust) the continentally derived insoluble particles, but only slightly more ¹⁰Be; uncertainties from curve-fitting typically are 20% of the ratios. Processes by which chemicals were transferred to the ice sheet were similar during the two periods, or at least had similar net effects (within on the order of 10%). Dry deposition was more important during the Younger Dryas than during the Preboreal hut wet deposition was more important than dry deposition for most species during cold as well as warm periods. Thus, contaminant concentration in the GISP2 ice core follow atmospheric concentrations more closely than do contaminant fluxes to the ice sheet. For the Younger Dryas to Preboreal climatic transition, changes in contaminant concentrations overestimate atmospheric changes by about 20% and changes in fluxes underestimate atmospheric changes by about 40% (Table 4).

Given sufficiently large data sets of chemical or particulate concentrations and annual-layer thicknesses in regions of simple ice flow, it should be possible to use the techniques outlined here to improve estimates of changes in atmospheric concentration of many chemical or particulate species in many locations and time intervals, and to provide time series with high time resolution as well as the low-resolution comparisons offered here. We anticipate that the reconstructed changes in atmospheric loadings will prove useful in a variety of paleoclimatic studies related to reconstructing wind speeds, source regions, changes in radiative scattering in the atmosphere, etc., and we are pursuing some such studies.