An oceanic origin for the increase of atmospheric radiocarbon during the Younger Dryas

1. Introduction

[2] Variation in Δ¹⁴C_atm is linked to climatic fluctuations, which affect the distribution of ¹⁴C throughout the Earth system, and to changes in geomagnetic and solar magnetic field intensity, which affects the production rate of ¹⁴C. While we are able to draw upon dendrochronology for Δ¹⁴C_atm after 12.4 kyr BP [Friedrich et al., 2004], prior to this we must rely on less well-constrained data mainly from marine archives, such as corals [Fairbanks et al., 2005] and foraminifera [Hughen et al., 2000]. Marine reconstructions assume a constant difference between the ¹⁴C age of the ocean and that of the contemporaneous atmosphere (the reservoir age, R) to calculate Δ¹⁴C_atm. A range of temporally invariant values of R, usually between 300 and 450 yr, have been assumed for marine Δ¹⁴C_atm reconstructions from different locations in the ocean. However, marine R can vary temporally and spatially as a result of atmospheric production rate or ocean ventilation changes [Butzin et al., 2005; Bondevik et al., 2006; Delaygue et al., 2003]. Indeed, comparison of tropical Atlantic data with a floating tree-ring record [Kromer et al., 2004] suggests a decrease in R of ∼200 years at the YD onset. Here, we use an Earth system model of intermediate complexity to examine the constant-R assumption, and reliability of ¹⁴C age calibration based on marine archives. In doing so, we explore potential causes of apparent variations in the marine radiocarbon record.

[3] We focus on a critical period that interrupted the last deglaciation: the Younger Dryas (YD; 12.9–11.7 kyr BP) cold reversal. According to Δ¹⁴C_atm reconstructions [Reimer et al., 2004; Fairbanks et al., 2005; Hughen et al., 2000], a large positive Δ¹⁴C_atm excursion occurred at the start of the YD (Figure 1a). Hughen et al. estimate an apparent increase of ∼50‰. The cause of this peak, and indeed the nature of Earth system change across the globe during this period are keenly debated [Goslar et al., 2000; Marchal et al., 2001]. This apparent excursion reached a peak Δ¹⁴C_atm value within 200 years of onset, and started concurrently with the initiation of cool Northern Hemisphere temperatures at the Bölling-Allerød/YD transition, as inferred from δ¹⁸O in ice cores [Rasmussen et al., 2006] (Figure 1b). Two key issues are: (1) the structure and amplitude of the Δ¹⁴C_atm anomaly in the absence of direct estimates of the atmospheric signal, and (2) the extent to which the Δ¹⁴C_atm anomaly can be attributed to changes in the production rate due to solar/geomagnetic variation and/or abrupt changes in Atlantic overturning, which influences atmosphere-ocean carbon exchange.

[4] The flux of the cosmogenic nuclide ¹⁰Be (Figure 1d) can be used as a proxy for production rate [Muscheler et al., 2004] because it is influenced by solar and geomagnetic intensity but, unlike ¹⁴C, is relatively unaffected by redistribution within the Earth system during changing climate conditions. Current ¹⁰Be records show no obvious excursion at the beginning of the YD, although uncertainties in ice accumulation rates [Muscheler et al., 2004] used to calculate ¹⁰Be production rate from ¹⁰Be concentration may mask variation. Alternatively, Pa/Th data from the North Atlantic provide evidence for reduced North Atlantic Deep Water (NADW) formation during the YD [McManus et al., 2004] which could affect Δ¹⁴C_atm via the associated reorganization of carbon reservoirs. Although Earth system models have previously simulated a Δ¹⁴C_atm increase in response to surface freshwater "hosing' and shutdown of NADW formation [Marchal et al., 2001; Butzin et al., 2005], they have significantly underestimated the apparent magnitude and rate of YD Δ¹⁴C_atm increase as derived from marine records.

[5] In this study, we examine the modelled spatial and temporal distribution of ¹⁴C within the various reservoirs of the global carbon cycle in response to both meltwater input and increased cosmogenic production using the GENIE-1 (Grid ENabled Integrated Earth system) intermediate-complexity model [Lenton et al., 2006]. The results were compared to empirical records to investigate the cause of the YD Δ¹⁴C_atm increase.

2. Model and Simulations

[6] GENIE-1 consists of a dynamic 3-D ocean, 2-D energy and moisture balance atmosphere, sea-ice, and terrestrial, ocean, and sediment carbon cycle components, on a 36 by 36 equal area grid with 16 vertical ocean levels [Lenton et al., 2006; Ridgwell et al., 2007; Ridgwell and Hargreaves, 2007]. It is forced with the annual average wind speed [Trenberth et al., 1989] and calibrated against present-day observations of surface air temperature and humidity, and 3D ocean temperature and salinity using a multi-objective tuning process [Price et al., 2006].

[7] For idealised Younger Dryas simulations, the model was prescribed with fixed ice sheet [Peltier, 1994] and orbital boundary conditions for 13 kyr BP. The model was first spun up for 40 kyr with a prescribed CO₂ concentration of 237 ppm and ¹⁴C production rate to produce Δ¹⁴C_atm of 210‰. A resulting AMOC maximum of 15.9 Sv and annual average SST of 14.8°C were obtained. During the actual Younger Dryas simulations, CO₂, ¹³C and ¹⁴C concentrations are predicted by the interactions between terrestrial, marine, and sediment carbon cycling and weathering.

[8] Firstly, the model was forced with an idealised freshwater flux in the high latitude Atlantic Ocean (Scenario A, Figure 2a), consistent with reconstructions of meltwater/iceberg discharge into the Atlantic and Arctic Ocean (Figure 1c) [Tarasov and Peltier, 2005] and equivalent to 7 m sea level rise in total, which is within the range of estimates from palaeo-records for the YD period [Fleming et al., 1998]. A second simulation was conducted in which an instantaneous increase in ¹⁴C production of 25% is prescribed for 400 yr (Scenario B, Figure 2a), which is the maximum magnitude of increase that could potentially be inferred from the ¹⁰Be record within uncertainties, to evaluate the response solely to production rate variation. In addition, a simulation was performed in which both a shutdown of NADW and idealised transient changes in production rate were prescribed, based on the increase in ¹⁴C production at the start of the YD derived from ¹⁰Be ice-core data (Scenario C, Figure 2a).

3. Model Results

[9] In scenario A, freshwater input is sufficient to shutdown NADW in the model, which takes ∼100 yr (Figure 2b). This produces an increase in Δ¹⁴C_atm of 30‰ over 800 yr (Figure 2c). The rate of Δ¹⁴C_atm increase is similar to that predicted in previous modelling studies [Marchal et al., 2001; Butzin et al., 2005]. However, this is significantly slower than the marine reconstructions apparently indicate (Figure 1a). Introducing a 25% increase in ¹⁴C production rate without freshwater-hosing (scenario B), generates a more rapid increase in Δ¹⁴C_atm (Figure 2c) more comparable with marine-reconstructed Δ¹⁴C_atm than the freshwater flux simulations. The final simulation (scenario C), which prescribed both freshwater-hosing and production rate changes, yielded Δ¹⁴C_atm increases of the same order as the other hosing experiments (∼40‰), but with a more rapid increase in Δ¹⁴C_atm than hosing alone.

[10] Reduction in ocean ventilation due to NADW shutdown (A) enables a longer period for exchange of ¹⁴C between atmosphere and Atlantic surface water. This results in surface water Δ¹⁴C values that are closer to the values in the atmosphere, and consequently lower R (Figure 3). Large reductions in R occur in the North Atlantic by year 200 (Figure 3e–3h). By model year 450 there is a considerable decrease in R over the whole of the Atlantic Ocean of up to 300 yr. Conversely, an increase in R is observed under regions such as the Arctic Ocean, where sea-ice extent is greater and, hence, air-sea gas exchange is reduced. The increase in R in the Arctic and North Atlantic is of the same order of magnitude as that seen in ocean core records [Bondevik et al., 2006]. Much smaller increases in R are modelled in the Pacific. This heterogeneous spatial distribution of R in scenario A is in contrast to that produced by an increase in ¹⁴C production rate, which produces a more uniform increase in R of up to 200 yr in both the Atlantic and Pacific (scenario B; not shown).

4. Model-Data Comparison

[11] Given the significant transient and spatial changes in modelled R, we next considered the implications for the robustness of the marine-based Δ¹⁴C_atm data reconstructions (which assume a constant R) at the YD onset. Using the original simulations A and B, we produced modelled time-series of ¹⁴C concentrations for the surface ocean in the same regions of the tropical Atlantic and Pacific as the paleo-archive (supplementary material, Figure S1) and use these to reconstruct Δ¹⁴C_atm by subtracting a value for R that is (wrongly) assumed to be constant [e.g., Delaygue et al., 2003]. When a constant R is imposed on the tropical Atlantic model surface ocean concentrations, the freshwater-hosing-only simulation results in an artificially large increase in reconstructed Δ¹⁴C_atm, which reaches a maximum value several hundred years before the model maximum Δ¹⁴C_atm (Figure 4a). When the constant-R reconstruction is performed on modelled Pacific surface ocean Δ¹⁴C, there is a significant delay before the rise in reconstructed Δ¹⁴C_atm begins (Figure 4a). In contrast to this, scenario B (¹⁴C production rate increase only), produces Δ¹⁴C_atm reconstructions that are damped in both the Pacific and Atlantic in comparison to the atmospheric values of Δ¹⁴C_atm (Figure 4b), and leads (Atlantic) and lags (Pacific) are not seen. In the production+hosing scenario (C), constant-R Δ¹⁴C_atm reconstructions (Figure 4c) produce a similar lead/lag pattern as Figure 4a, indicating the dominating effect of the change in ocean ventilation.

[12] To evaluate whether the phase and amplitude effects described above were present in the Δ¹⁴C_atm record, we compared smoothed marine-based Δ¹⁴C_atm records using data from the Atlantic [Fairbanks et al., 2005; Hughen et al., 2000] and Pacific [Fairbanks et al., 2005; Cutler et al., 2004; Burr et al., 2004] (Figure 4d). There is general agreement between the two ocean basins prior to the YD and during the latter half of the YD. however, at the start of the YD, the maximum Δ¹⁴C_atm occurs ∼300 yr earlier in the Atlantic reconstruction than the Pacific reconstruction, and is significantly larger. This pattern is comparable to the modelled scenario A (Figures 4a and 4c), and indicates that a decrease in Atlantic basin ventilation is the dominant contributing factor to the initial rise in reconstructed Δ¹⁴C_atm.

[13] Additional support for the model results is provided by a previous comparison of a floating tree ring record [Kromer et al., 2004] with the Cariaco Basin varved sediment Δ¹⁴C_atm record. Discrepancies between them suggested that reservoir ages in the Cariaco Basin (located in the Atlantic) decreased by ∼200 yr during the early part of the Younger Dryas, as found in our model simulations of the tropical Atlantic (Figure 3f).

5. Conclusions

[14] We demonstrate here that comparisons between marine ¹⁴C data from different oceans and Earth-system model simulations of Younger Dryas conditions can elucidate the mechanisms ¹⁴C variations. We find that the Younger Dryas Δ¹⁴C_atm increase is predominantly caused by ocean reorganisation rather than changes in production rate.

[15] For the first time, the possibility exists for reconciling empirical marine and modelled data for this period, and for understanding the causes of other variations in the radiocarbon record. Using our model-data comparison strategy, we conclude that marine-derived Δ¹⁴C_atm reconstructions that assume temporally invariant reservoir ages will yield erroneous phase relationships and amplitudes during other periods of rapid change in ocean overturning or ¹⁴C production that have occurred in the last 50 kyr. This has implications for the accuracy of the ¹⁴C age calibration efforts based on marine archives. However, there is potential for model simulations to inform and improve the error estimation of marine-based reconstructions of atmospheric radiocarbon abundance.