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 (Reference Davidson and LangwayDavidson, 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. Reference Davidson and LangwayDavidson, 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 10Be 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 (Reference Hodge, Wright, Bradley, Jacobel, Skou and VaughnHodge and others, 1990)). Modern mean annual temperature is about −31°C, modern mean accumulation is about 0.24 m year−1 of ice and melting occurs only about once per century or less frequently (Reference AlleyAlley and others, 1993; Reference Alley and AnandakrishnanAlley and Anandakrishnan, 1995).
Many of the data used here have been presented and discussed elsewhere (Reference AlleyAlley and others, 1993; Reference MayewskiMayewski and others, 1993a) The 10Be 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 Reference Mershon and ZielinskiMershon 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 −3), fluxes to the ice (molecules or particles m−2 year−1, calculated as the product of the concentration and the ice accumulation or flux in m year−1) and the age relative to the major late-glacial climatic events (e.g. Reference JohnsenJohnsen and others, 1992; Reference AlleyAlley and others, 1993; Reference TaylorTaylor 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 (Reference TaylorTaylor and others, 1992; Reference AlleyAlley 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 Reference Sehøtt, Waddington and RaymondSchøtt and others (1992) were exceedingly accurate, so we simply scaled the flow corrections from their model (Reference AlleyAlley and others, 1993). Further work on ice-flow corrections (Reference Bolzan, Waddington and AlleyBolzan and others, 1995; Reference Cutler, Raymond, Waddington, Meese and AlleyCutler and others, 1995) in general supports the Reference Sehøtt, Waddington and RaymondSchø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 (Reference TaylorTaylor and others, 1993a; Reference Alley, Gow, Johnsen, Kipfstuhl, Meese and ThorsteinssonAlley and others, 1995).
Data on concentrations of soluble major anions and cations (sulfate, nitrate, ammonium, chloride, sodium, potassium, magnesium and calcium) are from Reference MayewskiMayewski 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 (Reference Mershon and ZielinskiMershon 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.
10Be 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 10Be measured by accelerator mass spectrometry at the Lawrence Livermore National Laboratory Center for Accelerator Mass Spectrometry (Reference DavisDavis 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), 10Be (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 10Be 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 10Be 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 10Be, 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 10Be 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. Reference LegrandLegrand, 1987; Reference Davidson and LangwayDavidson, 1989; Reference Whung, Saltzman, Spencer, Mayewski and GundestrupWhung and others, 1994) relating contaminant flux to atmospheric concentration and ice flux is
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 1 is the dry-deposition velocity, k 2 is the dimensionless scavenging ratio for precipitation and Ki is the concentration of the contaminant species in the ice core.
Assume initially that k 1 and k 2 are constants (we discuss this assumption in the Model test section, below). Then Equation (1) shows that, if dry deposition dominates (k 2 Cib ≪ k 1 C i), the flux of some contaminant to the ice sheet will be proportional to its atmospheric concentration. If wet deposition dominates (k 2 Cib » k 1 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. Reference MayewskiMayewski and others, 1993a, Reference Mayewskib). 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 1 and k 2 are constants (see below), a plot of contaminant flux (ƒ i) vs ice accumulation (b) should produce a straight line with slope k 2 C i and intercept k 1 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, Ci , from its mean value, C i . Deviations may also be caused by measurement error, time variation of k 1 or k 2, 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 (Reference Raisbeck, Yiou, Jouzel, Petit, Barkov, Bard, Bard and BroeckerTaylor and others, 1992, Reference AlleyAlley 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 (Reference Whung, Saltzman, Spencer, Mayewski and GundestrupWhung 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 (Reference LegrandLegrand, 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. Reference TillTill, 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 1 and k 2 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, Riu: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
The two lines for the two time intervals have three rather than four free parameters: one slope, one intercept and the ratio Riu: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 Reference Press, Flannery, Teukolsky and VetterlingPress 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 10Be (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 10Be; Β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 10Be. 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 10Be 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 10Be 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 10Be 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 RcoulD: 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 RcoulD: 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 10Be data and show that any variations in k 1 and k 2 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 10Β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 10Be data.
10Be is a long-lived (1.5 × 106year) 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 (Reference Lal and PetersLal and Peters, 1967). Residence time is almost a year in the stratosphere but only about 3 weeks in the troposphere.
10Be 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 (Reference Beer, Andree, Oeschger and StaufferBeer and others, 1983, Reference Beer1990), events of increased production have occurred about 35 000 and perhaps 60 000 ABP, probably in response to supernovas or magnetic field events (Reference Raisbeck, Yiou, Bourles, Lorias, Jouzel and BarkovnRaisbeck and others, 1987, Reference Raisbeck, Yiou, Jouzel, Petit, Barkov, Bard, Bard and Broecker1992) and slow variations in production similar to those for 14C may have occurred in response to slow changes in the Earth’s magnetic field (Reference Lal, Bard and BroeckerLal, 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 (Reference McHargue and DamonMcHargue and Damon, 1991; Reference Lal, Bard and BroeckerLal, 1992; Reference Mazaud,, Laj, Bard, Arnold, Trie, Bard and BroeckerMazaud and others, 1992; Reference Raisbeck, Yiou, Jouzel, Petit, Barkov, Bard, Bard and BroeckerRaisbeck and others, 1992).
Figure 1 shows the time series of 10Be flux during the deglaciation and Figure 2a shows how 10Br 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 10Be flux in central Greenland.
Because of the likelihood that 10Be production has been nearly constant, we can use Figure 2 to test whether our model fits the data, and thus whether k 1 and k 2 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 10Be flux between cold (glacial or stadial) and warm (interglacial or interstadial) times must fall between rather narrow limits if k 1 and k 2 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 10Be is removed from the atmosphere by wet deposition (Reference McHargue and DamonMcHargue 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. Reference Kutzbach and GuetterKutzbach and Guetter, 1986; Reference Lautenschlager and HerterichLautenschlager and Herterich, 1990; Reference JoussaumeJoussaume, 1993). If 10 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 10Be 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 1 and k 2, 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 10Be is not long enough for complete global mixing in the troposphere. Suppose instead that no lateral mixing of 10Be occurs in the atmosphere, so that the flux of 10Be 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 10Be 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 1 and k 2 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 (Reference Lal and PetersLal 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. Reference Rind, Peteet, Broecker, Mclntyre and RuddimanRind and others, 1986). Thus, we expect that 1 < R coulD: warm 10Be < 1.5 for 10Be between the Younger Dryas and Preboreal times, if k 1 and k 2 are constants.
Thus far, our data set includes many 10Be 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. Reference TillTill, 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 10Be 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 1 and k 2 for 10Be 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 10Be 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 10Be 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 Reference McHargue and DamonMcHargue 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 10Be-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. Reference ThompsonThompson, 1977; Reference Petit, Briat and RoyerPetit and others, 1981; Reference Mershon and ZielinskiMenhon 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. Reference WanneckWarneck, 1988; Reference SchumannSchumann, 1991; Reference Hillamo, Kerminen, Maenhaut, Jaffrezo, Balachandran and DavidsonHillamo 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. Reference Bergin, Jaffrezo, Davidson, Caldow and DibbBergin and others, 1994; Reference Dibb, Talbot and BerginDibb 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 1 and k 2), 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. Reference Hillamo, Kerminen, Maenhaut, Jaffrezo, Balachandran and DavidsonHillamo 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 Reference WanneckWarneck, 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. Reference WanneckWarneck, 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 (Reference Cunningham and WaddingtonCunningham and Waddington, 1993; Reference Hillamo, Kerminen, Maenhaut, Jaffrezo, Balachandran and DavidsonHillamo 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 Reference WanneckWarneck (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. Reference Clausen, Langway, Oeschger and LangwayClausen and Langway, 1989; Reference Delmas, Legrand, Oeschger and LangwayDelmas and Legrand, 1989; Reference Mayewski, Spencer, Twickler and WhitlowMayewski 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 Reference Clausen, Langway, Oeschger and LangwayClausen and Langway, 1989; Reference Delmas, Legrand, Oeschger and LangwayDelmas and Legrand, 1989; Reference Mayewski, Spencer, Twickler and WhitlowMayewski 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 (Reference Clausen, Langway, Oeschger and LangwayClausen and Langway, 1989; Reference Delmas, Legrand, Oeschger and LangwayDelmas and Legrand, 1989; Reference Mayewski, Spencer, Twickler and WhitlowMayewski and others, 1990; Reference Legrand, Angelis, Staffelbach, Neftel and StaufferLegrand 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 (Reference LegrandLegrand, 1987: Reference Davidson and LangwayDavidson, 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, Reference Davidson, Santhanam, Fortmann and OlsonDavidson and others, 1985; Reference Dibb, Jaffrezo and LegrandDibb and others, 1992; Reference Jaffrezo, Hillamo, Davidson and MaenhautJaffrezo 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. Reference Peel, Mulvaney and DavisonPeel and others, 1988; Reference Cuffey, Alley, Grootes and AnandakrishnanCuffey and others, 1992, Reference Cuffey, Alley, Grootes, Bolzan and Anandakrishnan1994).
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, Reference LegrandLegrand (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 Reference Peel, Mulvaney and DavisonPeel 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 (Reference Bergin, Jaffrezo, Davidson, Caldow and DibbBergin 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. Reference TillTill, 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 Reference MayewskiMayewski 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 10Be data indicates that model errors are small for 10Be 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 10Be but not for other species. We note that we chose to test the model against 10Be, and we estimated the limiting behaviour of 10Be, before we had analyzed the 10Be data; thus, any claim for uniqueness of 10Be 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. Reference Petit, Briat and RoyerPetit and others, 1981; Reference Herron, Langway, Langway, Oeschger and DansgaardHerron and Langway, 1985; Reference HarveyHarvey, 1988; Reference Thompson, Davis, Mosley-Thpmpson and LiuThompson and others, 1988; Reference MayewskiMayewski 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. Reference Davidson and LangwayDavidson, 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 10Be; 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.
Acknowledgements
We thank the GISP2 Science Management Office, the U.S. 109th Air National Guard, the Polar Ice Coring Office and the U.S. National Ice Core Laboratory for scientific and logistical support, the U.S. National Science Foundation Office of Polar Programs for financial support (R.B.A. also thanks NASA-EOS and the D. and L. Packard Foundation for financial support), and R. Bales, M. Bergin, J. Dibb, D. Lamb, D. Meeker, E. Saltaman, E. Wolff an anonymous reviewer, and other GISP2 colleagues for advice and helpful suggestions. This work was partially supported under the auspices of the U. S. Department of Energy by Lawrence Livermore National Laboratory under contract W-7405-Eng-48.
