A South Atlantic island record uncovers shifts in westerlies and hydroclimate during the last glacial

Changes in the latitudinal position and strength of the Southern Hemisphere westerlies (SHW) are thought to be tightly coupled to important climate processes, such as cross-equatorial heat fluxes, Atlantic Meridional Overturning Circulation (AMOC), the bipolar seesaw, Southern Ocean ventilation and atmospheric CO₂ levels. However, many uncertainties regarding magnitude, direction, and causes and effects of past SHW shifts still exist due to lack of suitable sites and scarcity of information on SHW dynamics, especially from the last glacial. Here we present a detailed hydroclimate multiproxy record from a 36.4–18.6 kyr old lake sediment sequence on Nightingale Island (NI). It is strategically located at 37^∘ S in the central South Atlantic (SA) within the SHW belt and situated just north of the marine Subtropical Front (SF). This has enabled us to assess hydroclimate changes and their link to the regional climate development as well as to large-scale climate events in polar ice cores. The NI record exhibits a continuous impact of the SHW, recording shifts in both position and strength, and between 36 and 31 ka the westerlies show high latitudinal and strength-wise variability possibly linked to the bipolar seesaw. This was followed by 4 kyr of slightly falling temperatures, decreasing humidity and fairly southerly westerlies. After 27 ka temperatures decreased 3–4 ^∘C, marking the largest hydroclimate change with drier conditions and a variable SHW position. We note that periods with more intense and southerly-positioned SHW seem to be related to periods of increased CO₂ outgassing from the ocean, while changes in the cross-equatorial gradient during large northern temperature changes appear as the driving mechanism for the SHW shifts. Together with coeval shifts of the South Pacific westerlies, our results show that most of the Southern Hemisphere experienced simultaneous atmospheric circulation changes during the latter part of the last glacial. Finally we can conclude that multiproxy lake records from oceanic islands have the potential to record atmospheric variability coupled to large-scale climate shifts over vast oceanic areas.

1 Introduction

The SHW are a major determinant of hydroclimate in the Southern Hemisphere (SH). In coupling marine and atmospheric processes, they are thought to have played a pivotal and multifaceted role during and at the end of the last ice age by triggering changes in ocean–atmosphere CO₂ fluxes through physical processes (Saunders et al., 2018; Toggweiler and Lea, 2010) and Fe fertilization of the Southern Ocean through varying dust deposition (Lamy et al., 2014; Martin and Fitzwater, 1988; Martínez-García et al., 2014), as well as regulating the salt and heat leakage from the Agulhas Current to the Atlantic Meridional Overturning Circulation (AMOC; Bard and Rickaby, 2009). In addition, changes in AMOC, SHW strength and position, and Southern Ocean upwelling seem to have been important mechanisms for different glacial CO₂ modes (Ahn and Brook, 2014). The position of the SHW during glacial times is debated with some arguing for a northwards displacement (Toggweiler et al., 2006) while others argue for a southwards move (Sime et al., 2013, 2016) during the Last Glacial Maximum (LGM), relative to the present. Holocene data also suggest an expanding–contracting SHW zone (Lamy et al., 2010). With these multiple scenarios, the pattern of SHW shifts and their detailed role for ocean ventilation and the global carbon cycle remains unclear. It is postulated that the SHW moved in concert with rapid climate shifts recorded in Greenland ice cores known as Dansgaard–Oeschger (DO) cycles (Markle et al., 2016) and that these shifts are part of interhemispheric climate swings involving heat exchange between the hemispheres through the atmosphere and the ocean, with atmospheric heat fluxes partly compensating anomalous marine heat fluxes (Pedro et al., 2016). Whether SHW zonal shifts only occurred in the Pacific sector of the Southern Ocean (Chiang et al., 2014) or if they occurred throughout the SH is another crucial question (Ceppi et al., 2013). Other key climate issues relate to the effects and areal extent of the bipolar seesaw mechanism (Broecker, 1998; Stocker and Johnsen, 2003) and any signs of an early and long temperature minimum at southern mid-latitudes matching Antarctic LGM (EPICA Community Members et al., 2006). The lack of climate proxy records directly reflecting atmospheric conditions in the central South Atlantic means that such information at these latitudes during the glacial are primarily based on remote proxy records or climate model simulations. This results in a largely unconstrained understanding of glacial conditions over vast parts of the mid-South Atlantic, especially between 20 and 50^∘ S where archives reflecting atmospheric processes are absent. In this study our aim is to use hydroclimate and temperature proxy records as well as climate model output to reconstruct changes in the position of the SHW in the Atlantic sector, reconstruct hydroclimate changes, and identify interhemispheric linkages including the bipolar seesaw during the last glacial period. We also explore links between past SHW strength and atmospheric CO₂. For these questions the Tristan da Cunha archipelago is uniquely situated in the South Atlantic (Fig. 1a).

Figure 1 (a) The position of the Tristan da Cunha island group in the South Atlantic and the 1000 mbar mean annual wind speed (m s⁻¹) for 1980–2010 according to NCEP/NCAR reanalysis data showing yellow to red colors for the zone of the Southern Hemisphere westerlies, and the positions of the Subtropical Front (SF) and the Polar Front (PF) as dashed lines. (b) The three main islands of the Tristan da Cunha island group. (c) The position and size of the four overgrown lake basins, so-called ponds (labeled 1P–4P), on Nightingale Island with 100 m contour lines.

2 Study site

The Tristan da Cunha island group (TdC) at 37.1^∘ S (Fig. 1) sits strategically at the northern boundary of the SHW (Fig. 1a), a few degrees north of the Subtropical Front (SF), where sea surface temperatures (SSTs) and salinities decrease by 3–4 ^∘C and 0.3 ‰, respectively. Annual mean air temperature and precipitation are 14.3 ^∘C and approximately 1500 mm, respectively, with highest precipitation in austral winter when the SHW impact is largest. The record presented here is from the 1st Pond (labeled 1P), an overgrown crater lake (200 m×70 m, 207 m a.s.l.) today forming a peat bog in the central part of Nightingale Island (NI) (Figs. 1 and 2), a volcanic island dominated by trachytic bedrock. Its drainage area is about twice the size of the peat bog and is thus sensitive to changes in the precipitation and evaporation balance (P/E). Previous studies from NI show that the area experienced shifts in precipitation during the Holocene (Ljung and Björck, 2007) and partly also during the Last Glacial Termination (Ljung et al., 2015), mainly attributed to the changing impact of the SHW. These data also indicate a southerly displacement of the Intertropical Convergence Zone (ITCZ) during the Heinrich 1 event (H1) and warming in the South Atlantic as a consequence of reduced AMOC, causing the lake basin to dry out, and create a hiatus between 18.6 and 16.2 ka (Ljung et al., 2015). Here we present a multiproxy study of the sediments that accumulated before this hiatus dating to 36.4–18.6 ka, covering the younger part of Marine Isotope Stage 3 (MIS 3) and most of MIS 2, a climatically very dynamic period with Antarctic Isotope Maxima (AIM), DO events and H events. In spite of its fairly northern position in relation to Antarctica, we hypothesize that TdC was impacted by such events in terms of shifts of SHW, which we aim to test by using a suite of proxies.

Figure 2 Photograph from Nightingale Island. The overgrown lake basins of 1st Pond and 2nd Pond are shown, with the higher-situated 1st Pond in the background, seen towards the southeast. Note the albatross chicks (white dots) and the four people on 2nd Pond for scale. Photo Svante Björck.

3 Material and methods

A large set of proxy data were analyzed, including chemical (N, XRF (X-ray fluorescence) elemental concentrations and isotopes (¹³C, ¹⁵N, ²H or D)), biological (TOC (total organic carbon), molecular fossils such as n-alkanes, glycerol dialkyl glycerol tetraether lipids (GDGTs), pollen and diatom assemblages, and biogenic silica (BSi)), and physical (magnetic susceptibility (MS)) parameters. Some proxies provide information about local changes such as soil conditions and erosion (C∕N ratios, ¹³C and MS), weathering (major element data), vegetation composition (pollen, n-alkane distributions), organic productivity (TOC and BSi), lake conditions and levels (diatoms, BSi, δD values of short-chained n-alkanes), and bird impact (¹⁵N). Others display regional changes in hydroclimate, such as mean annual air temperature (MAAT) and mean summer air temperature (MST) from the GDGT lipids and the source water of terrestrial and aquatic plants including evaporative conditions (hydrogen isotopes, δD). Observations of the isotopic content of precipitation are very sparse around TdC, and therefore we have investigated the hydroclimate variability with an isotope-enabled climate model. In addition, we have performed principal component analysis (PCA) to distinguish the influence of the different proxies on samples (see Methods, Sect. 3.13). Most of our data are found in the Supplement.

3.1 Field work, handling of cores and sample collection

Two weeks of field work on NI were carried out in February 2010 and drilling was carried out using Russian chamber samplers providing 1 m long cores (∅=50 and 75 mm) with overlaps of 15–50 cm between each cored section. The ketch Ocean Tramp provided the transport from the Falkland Islands to TdC and back to Uruguay. In order to penetrate as deep as possible into the very stiff sediments, a chain hoist was used for coring the deeper parts of the sequences. The sediments were described immediately in the field before being wrapped in plastic film and PVC tubes. Upon arrival in Uruguay, the cores were transported to the Geology Department in Lund where they were stored in a cold room. Before subsampling for the different proxy analyses, the field-based lithostratigraphy and correlations between individual core sections were adjusted in the laboratory. This was aided by magnetic susceptibility (κ) measurements, which give a relative estimate of the magnetic mineral concentration, to confirm and adjust the visual correlation between overlapping core segments.

3.2 Radiocarbon dating and age model

The radiocarbon-dated material consisted of 1 cm thick, organic-rich, bulk sediment. All 41 dated samples were pre-treated and measured at the Lund University Radiocarbon Dating Laboratory with single-stage accelerator mass spectrometry (SSAMS). The age model (Fig. 3) was constructed using the OxCal software package (Bronk Ramsey, 1995, 2009a). To minimize subjective user input, we ran the age model with a general outlier model (Bronk Ramsey, 2009b) and a variable k value that lets the model itself determine the sedimentation rate variability (Bronk Ramsey, 2008). For calibration we use the Southern Hemisphere calibration data set, SHCal13 (Hogg et al., 2013).

3.3 Measurements for magnetic susceptibility

Magnetic susceptibility (κ) was measured using a Bartington MS2E1 high-resolution surface scanning sensor coupled to a Tamiscan automatic logging conveyor. Measurements were carried out on non-sampled half cores and with a resolution of 5 mm and with results shown in 10⁻⁶ SI units. The magnetic susceptibility gives a relative estimate of the ability of the material to be magnetized, i.e., the magnetic mineral concentration.

Figure 3 Age model for the sediments at 1st Pond, Nightingale Island. Top plot: radiocarbon-based age–depth model (black lines encompass the 68.2 % probability interval). The patches indicate the calibrated probability distributions of each radiocarbon date for un-modeled (single) dates (gray patch) and their posterior distributions when modeled as a P sequence: green patches indicate agreement indices of > 60 % and red patches agreement indices of < 60 %, i.e., outliers. Bottom plot: sediment accumulation rates (SAR, mm a⁻¹) based on the mean age–depth model shown in the top panel.

3.4 XRF analyses

A handheld Thermo Scientific portable XRF analyzer (h-XRF) Niton XL3t 970 GOLDD+ set to the Cu/Zn mining calibration mode was used. The instrumentation provides highly accurate determinations for major elements (Helfert et al., 2011). All analyses were performed on freeze-dried sediments from the 1P cores using an 8 mm radius spot size in order to obtain representative values. The elemental detection depends partly on the duration of the analysis at each point; this is especially true for the lighter elements such as Mg, Al, Si, P, S, Cl, K and Ca. For this reason the measurement time of each sample was set to 6 min. Although a larger suite of elements was acquired, we have chosen to work with Al, Si, P, S, K, Ca, Ti, Mn, Fe, Rb, Sr and Zr. These elements were selected based on their analytical quality (i.e., level above the detection limit) and with the help of principal component analysis (PCA). PCA was made using JMP 10.0.0 software in correlation mode using a Varimax rotation. Before analysis all data were converted to Z scores calculated as $(X_i-X_{\mathrm{avg}})/X_{\mathrm{std}}$, where X_i is the normalized elemental peak areas and X_avg and X_std are the series average and standard deviation, respectively, of the variable X_i. A Varimax rotation allocates the component variables which are highly correlated (sharing a large proportion of their variance) – imposing some constrains in defining the eigenvectors. By grouping together elements showing similar variation, the chemical signals tend to be clearer and key elements are better identified. To simplify the interpretation of our principal components (PCs) we employ a modified Chemical Index of Alteration (CIA), see Fig. 4d, as defined by Nesbitt and Young (1982): CIA $=\left[{\mathrm{Al}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{3}}/({\mathrm{Al}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{3}}+\mathrm{CaO}+\mathrm{NaO}+{\mathrm{K}}_{\mathrm{2}}\mathrm{O})\right]\times \mathrm{100}$. This index expresses the relative proportion of Al₂O₃ to the more labile oxides and is an expression of the degradation of feldspars to clay minerals. Since we have no NaO data, we call it a modified CIA.

3.5 C and N analyses

Dried and homogenized samples every 1–2 cm were analyzed with a Costech Instruments ECS 4010 elemental analyzer. The accuracy of the measurements is better than ±5 % of the reported values based on replicated standard samples. To account for C∕N atomic ratios, the ratio was multiplied by 1.167.

3.6 ¹³C and ¹⁵N analyses

Dried, homogenized bulk samples were measured using a Thermo Fisher DeltaV ion ratio mass spectrometer. The isotopic composition of samples is reported as conventional δ values in parts per thousand relative to the Vienna Pee Dee Belemnite (¹³C) and atmospheric ¹⁴N (¹⁵N): δ_sample (‰) = [(R_sample−R_standard)/(R_standard)] ×1000, where R is the abundance ratio of ¹³C∕¹²C in the sample or in the standard.

3.7 Pollen analyses

Sixty-four levels were subsampled and analyzed for their pollen content. Pollen samples of 1 cm³ were processed following standard method A as described by Berglund and Ralska-Jasiewiczowa (1986) with added Lycopodium spores for determination of pollen concentration values. Counting was made under a light microscope at magnifications of 400× and 1000×. The aim was to count at least 500 pollen grains in every sample, which was almost achieved (mean sum of 565 pollen grains and mean sum of 870 pollen grains and spores). Identification of pollen grains and spores was facilitated by published photos (Hafsten, 1960), standard pollen keys (Moore et al., 1991) and a small collection of type slides from Tristan da Cunha borrowed from the University Museum of Bergen, the Natural History Collections. The pollen percentage diagram (Fig. 5) was plotted in C2 software (Juggins, 2007). Warm/cold pollen ratios were calculated as ($W_{\mathrm{p}}/W_{\mathrm{p}}+C_{\mathrm{p}}$), where warm pollen types (W_p) are from plants only found below 500 m a.s.l. and cold pollen types (C_p) are from plants only found above 500 m a.s.l.

3.8 Diatom analyses and diatom environmental ratios

One hundred and seventy-nine levels of 0.5 cm thick sediment segments were subsampled to analyze their diatom content. For preparation of diatom slides, ∼200 mg of freeze-dried sediment was oxidized with 15 % H₂O₂ for 24 h, then treated with 30 % H₂O₂ for a minimum of 24 h, and finally heated at 90 ^∘C for several hours. A known quantity of DVB (divinylbenzene) microspheres was added to 200 µL aliquots of the digested and cleaned slurries in order to estimate diatom concentrations (Battarbee and Keen, 1982). The diatoms were mounted in Naphrax^® medium (refractive index =1.65). Three hundred valves or more per sample were counted in most samples and identified largely using published diatom floras (Krammer and Lange-Bertalot, 1986; Le Cohu and Maillard, 1983; Moser et al., 1995; Van de Vijver et al., 2002). Diatom results are expressed as relative percent (%) abundance of each taxon (Fig. S3 in the Supplement) and also as total concentrations of valves per gram of dry sediment.

Freshwater diatom species are excellent indicators of water quality, particularly of pH, conductivity and dissolved nutrients (Battarbee et al., 2001). Sedimentary diatom assemblages inter alia can be used to reconstruct past changes in water quality using the ecological indicator information for each species. Where suitable modern diatom–water quality calibration data sets exist, transfer functions can be generated to reconstruct these changes. However, in sediment records where diatom diversity is low and affinities of some species are not firmly established, placing diatom taxa into ecological or environmental preference groups using literature attributions and field experience can be used to generate ratio scores relevant to past conditions. The 1st Pond assemblages are suitable for such an approach, particularly for inferring changes in habitat and water acidity. The acid diatom index ratio is derived from the sum of acid water indicating taxa comprising Aulacoseira, Frustulia, Pinnularia and Eunotia compared to that of the fragilarioid tychoplanktonic taxa. Proportions of acidity-tolerant to acidity-intolerant diatom taxa indicate water pH; total tychoplankton (temporary phytoplankton) vs. total benthic taxa relate to open-water conditions; subaerial and terrestrial taxa vs. the total assemblage indicate wetland development and/or in-washed material.

3.9 Biogenic silica analyses

The 310 samples were analyzed using a wet-alkaline digestion technique (Conley and Schelske, 2001). Samples were freeze-dried and gently ground prior to analysis. Approximately 30 mg of sample was digested in 40 mL of a weak base (0.47 M Na₂CO₃) at 85 ^∘C for a total duration of 3 h. Subsamples of 1 mL were removed after 3 h and neutralized with 9 mL of 0.021 M HCl. Dissolved Si concentrations were measured with a continuous-flow analyzer applying the automated molybdate blue method (Grasshoff et al., 1983). Biogenic silica content in lake sediments is a rough proxy for lake productivity.

3.10 Lipid biomarker and compound-specific hydrogen isotopic analyses

The hydrogen isotopic composition (δ notation) of n-alkanes was analyzed by gas chromatography isotope ratio monitoring mass spectrometry (GC-IRMS) using a Thermo Finnigan Delta V mass spectrometer interfaced with a Thermo Trace GC 2000 using a GC Isolink II and Conflo IV system. Helium was used as a carrier gas in constant flow mode and the compounds separated on a Zebron ZB-5HT Inferno GC column (30 m × 0.25 mm × 0.25 µm). Lipid extraction was performed on freeze-dried samples by sonication with a mixture of dichloromethane and methanol (DCM–MeOH 9:1 v∕v) for 20 min and subsequent centrifugation. The process was repeated three times and supernatants were combined. Aliphatic hydrocarbon fractions were isolated from the total lipid extract using silica gel columns (5 % deactivated) that were first eluted with pure hexane (F1) and subsequently with a mixture of DCM–MeOH (1:1 v∕v) to obtain a polar fraction (F2). A saturated hydrocarbon fraction was obtained by eluting the F1 fraction through 10 % AgNO₃−SiO₂ silica gel using pure hexane as eluent. The saturated hydrocarbon fractions were analyzed by gas chromatography mass spectrometry for identification and quantification, using a Shimadzu GCMS-QP2010 Ultra. C₂₁ to C₃₃ n-alkanes were identified based on mass spectra from the literature and retention times. The concentrations of individual compounds were determined using a calibration curve made using mixtures of C₂₁–C₄₀ alkanes of known concentration. More details about the GC-IRMS method, including GC oven temperature program, instrument performance and reference gases used, are given in Yamoah et al. (2016). The average standard deviation for δD values was 5 ‰. Due to low sea levels during the time period of our proxies, the δD values of the n-alkanes were ice volume corrected (Tierney and deMenocal, 2013), δD${}_{\mathrm{corr}}=({\mathit{\delta }\mathrm{D}}_{\mathrm{wax}}+\mathrm{1000})/({\mathit{\delta }\mathrm{O}}_{\mathrm{w}}^{\mathrm{18}}\times \mathrm{8}\times \mathrm{0.001}+\mathrm{1})-\mathrm{1000}$, with interpolated ocean water δO${}_{\mathrm{w}}^{\mathrm{18}}$ values (Waelbroeck et al., 2002).

Isoprenoid and branched glycerol dialkyl glycerol tetraethers (GDGTs) were measured on the F2 fractions after filtration through 0.45 µm PTFE filters and reconstitution into a known volume of methanol. Analysis was done using a Thermo–Dionex HPLC (high-performance liquid chromatography) connected to a Thermo Scientific TSQ Quantum Access triple quadrupole mass spectrometer, using an APCI (atmospheric pressure chemical ionization) interface. Chromatographic separation was achieved using a reverse phase method similar to the one used by Zhu et al. (2013). Partially coeluting GDGT isomers were integrated as one peak in order to obtain data comparable to the normal phase method that has been in use by the community since Weijers et al. (2007).

A basic prerequisite for the valid use of brGDGTs (branched GDGTs) is a relatively high branched-over-isoprenoid tetraether (BIT) index, which was 1.00 throughout the core. Reconstructed pH values, based on the cyclization of branched tetraether (CBT) ratio (Weijers et al., 2007), were stable at 6.6±0.1 over the length of the core, which means that temperature is the dominant environmental factor exerted on the brGDGT distribution. At the time of measurement, we had not adopted the new method which separates between 5-methyl and 6-methyl branched GDGTs (De Jonge et al., 2014). As a consequence, we do not have individual quantifications of 5-methyl and 6-methyl branched GDGT isomers needed to use the revised 5-methyl branched tetraether (MBT_5me) temperature proxy for mineral soils (De Jonge et al., 2014) or peat (Naafs et al., 2017), which gives lower root-mean-square error than the original terrestrial (soil) calibration (Weijers et al., 2007). However, since our data are from lake sediments, we argue that GDGT-based temperature proxy calibrations based on lake surveys is in any case a more valid approach. Indeed, using the original temperature calibration of Weijers et al. (2007) based on soils resulted in very low temperatures between 0 and 6 ^∘C, a cold bias observed in other studies from lakes. This bias is probably due to the addition of in situ produced brGDGTs on top of any brGDGTs eroded from land (Loomis et al., 2012; Pearson et al., 2011). Our record could be biased by a changing ratio of soil- and lake-derived GDGTs, where a greater relative contribution of terrestrial-derived GDGTs would result in a warm bias if a lake calibration is used. However, we do not find a correlation between GDGT-derived temperature and two proxies for terrestrial influx, the C∕N ratio and magnetic susceptibility, but rather the opposite. We used two lake calibration sets: (a) the one of Pearson et al. (2011), based on a global lacustrine data set and using mean summer temperatures (MST), including samples from nearby South Georgia island in the South Atlantic, and (b) a calibration based on a large data set of east African lakes from different altitudes (Loomis et al., 2012), using mean annual air temperatures (MAAT), and which is also applicable outside of east Africa (Loomis et al., 2012). It is impossible to test which of these two proxy records would reflect past conditions more accurately. However, the two reconstructions strongly covary, with a difference between reconstructed MST and MAAT of approximately 5 ^∘C.

3.11 Calculation of insolation values

A long-term numerical solution for earth's insolation quantities (Laskar et al., 2004) was used for the insolation values, 37–18 kyr at 37^∘ S, and calculated with the AnalySeries program. While the austral winter values were based on mean daily June–August insolation (W m⁻²), the mean austral summer values were based on the mean daily December–February insolation.

3.12 Isotope model simulation

The isotope model analysis is based on a 1200-year simulation using the isotope-enabled version of the ECHAM5/MPI-OM earth system model (Werner et al., 2016) run with natural and anthropogenic forcings for 800 to 2000 CE (Sjolte et al., 2018). Horizontal resolution of the atmosphere is $\mathrm{3.75}{}^{\circ }\times \mathrm{3.75}{}^{\circ }$ (T31) with 19 vertical layers, while the ocean has a horizontal resolution of $\mathrm{3}{}^{\circ }\times \mathrm{1.8}{}^{\circ }$ with 40 vertical layers. The model includes isotope fractionation for all phase changes in the hydrological cycle, including below-cloud evaporation. Since both the present-day situation and our Nightingale Island record show a continuous impact from the westerlies, we deem it valid to use this late Holocene simulation as an analog for interpreting the variability in the westerlies during the time period of study. The outcome of the simulation is presented in the results section but further investigation of the model run shows that the multidecadal variability in δD at TdC is related to the phase of the Southern Annular Mode, indicating that isotopic variability at TdC is sensitive to large-scale SH climate variability (Fig. S4).

Figure 4 Antarctic ice core data and some proxy data from the sediments in 1st Pond between 36.4 and 18.6 ka. (a) The EPICA Dome Dronning Maud Land (EDML) δ¹⁸O record (EPICA Community Members, 2006) showing AIM7–2. (b) GDGT-based mean annual air temperature (MAAT) and mean summer temperature (MST) calibrated to Pearson et al. (2011) and Loomis et al. (2012), respectively. (c) Warm pollen ratios, percent (%) Phylica arborea pollen and percent (%) Ophioglossum spores. (d) Modified chemical index of alteration (CIA). (e) Percent (%) Cyperaceae pollen and n-C₂₃ alkanes. (f) δD values (‰) of n-C₂₃ and n-C₂₇₋₃₁ alkanes. (g) Acid diatom ratios. (h) Percent (%) terrestrial diatoms. (i) Open water diatom ratios. (j) Percent (%) biogenic silica (BSi). (k) C∕N ratios. (l) % total organic carbon (TOC). (m) Magnetic susceptibility (MS) expressed as 10⁻⁶ SI units. All proxies relate to the age scale on the x axes. Note the two thick gray lines (31 and 26.5 ka) indicating the position of the three PCA zones (Fig. 6).

3.13 Principal component analysis (PCA)

PCA was performed with 14 of our proxies (Fig. 6b) that we expect to respond to hydroclimate changes but without the MAAT values in order to test the MAAT values vs. other climate proxies, using the C2 program (Juggins, 2007). The aim was to display the impact of different combinations of proxies on the samples in a biplot (Fig. 6b), as discussed below in Sect. 4.2. All proxy data were centered and standardized before calculation.

Figure 5 Pollen diagram from 1st Pond, Nightingale Island. The diagram shows relative abundance (%) of the pollen taxa. Note that it is both related to depth (cm) and age (ka) on the y axis, the latter according to the age–depth model in Fig. 3.

Figure 6 Principal component analysis (PCA) of 14 proxies from 1st Pond. (a) Scores of the first two principal components related to age and the three PCA zones. Note that negative values point upwards and how PC1 and PC2 values are related to temperatures and SHW to the right y axis. (b) PCA plot shows the loadings of the 14 proxies (shown as red dots and black text with reference to proxies in Fig. 4, except for Empetrum rubrum). PC1 (red-brown) and PC2 (blue) account for 38.1 % and 13.4 % of the variance, respectively. The interpretations of the two axes are shown by red-brown and blue texts, and the interpretations of the four segments are based on the combined positions of the proxies in the plot (shown in four different colors).

4 Results

4.1 An island record of glacial climate in the central South Atlantic

Thirty-nine 1 m long overlapping cores were taken in February 2010 from three overgrown crater lakes (Fig. 1c) between lava ridges (Anker Björk et al., 2011). The 1st Pond (1P) was exceptional in that it was the only site where sediments older than 18.6 kyr were recovered. At 1P the 16.2–18.6 kyr hiatus (Ljung et al., 2015) is marked by a thin silt lamina at 618.8 cm. We retrieved five overlapping cores below the hiatus with 318.2 cm of sediments before coring was obstructed at 937 cm by suspected bedrock or boulders. These cores were correlated by lithology and magnetic susceptibility (MS). The lower 162 cm consists of a dark-brown slightly silty gyttja, overlain by a gray-brown silty clay gyttja, all deposited under anaerobic conditions. Because of the low concentration of plant macro-fossil remains our chronology is based on 41 ¹⁴C dates of 1 cm thick bulk sediment samples between 620 and 936 cm (Table S1 in the Supplement). Comparisons of ¹⁴C dates of bulk sediment and plant remains (wood and peat) have shown good concordance (Ljung et al., 2015; Ljung and Björck, 2007), and the most likely explanation for the seven clear outliers (Fig. 3) is possibly a combination of statistical noise and contamination from small amounts of recent material. Our age model displays a mean sedimentation rate of 0.18 mm yr⁻¹ but with considerable variation.

Figure 7 Zonal mean changes in wind, temperature and precipitation related to δD variability at TdC. (a) Time series of simulated precipitation-weighted annual mean δD at TdC, with values above and below the 95th percentile indicated with red and blue circles, respectively. This selection of high and low δD is used to define the data in (e–g). (b–d) Annual modeled (black) South Atlantic zonal mean (30 to 0^∘ W) westerly wind speed (850 mbar U wind, positive towards east), 2 m temperature (T2m) and precipitation compared to the 20th Century Reanalysis climatology 1981–2010 (gray) (Compo et al., 2011). (e–g) Composites of annual modeled zonal mean (30 to 0^∘ W) westerly wind speed (850 mbar U wind), 2 m temperature (T2m) and precipitation for high (red) and low (blue) δD at TdC defined in (a). (h–j) High minus low anomalies of model output are shown in (e)–(f). Circles indicate significant anomalies (p < 0.01) calculated using two-tailed Student's t test. The vertical bars in (b)–(j) show the latitude of NI at 37^∘ S.

In agreement with the supposed minimum age of pond formation through volcanic activity (Anker Björk et al., 2011), the bottom of 1P has an age of 36.4±0.3 ka. Our temperature records (Fig. 4b) show an oscillating pattern, with the largest change at 27.5 ka, and share similarities with the low-frequency variability in the EPICA Dome Dronning Maud Land (EDML) curve (Fig. 4a). Before 27.5 ka, MAAT and MST vary between 17 and 12 ^∘C and between 21 and 17 ^∘C, respectively, while the variation is between 13 and 9 ^∘C and between 18.5 and 15.5 ^∘C, respectively, after 27.5 ka. In terms of pollen as a local temperature indicator, it is known that Phylica arborea, Acaena sarmentosa and two Asteraceae plant types are sensitive to cold conditions (Ryan, 2007). They make up warm pollen types at NI, and the warm/cold pollen-types ratio (Fig. 4c) shows large variations until 31.4 ka, followed by a two-step decline (at 31.2 and 26.5 ka) largely in contrast to the spore abundance of the cold tolerant Ophioglossum opacum fern, and with a trend similar to the temperature curves. In comparison to Holocene sediments from NI (Ljung and Björck, 2007), the glacial pollen record from 1P (Fig. 5) shows less variability, and the most distinct difference is the very low abundance of the only tree species pollen on the island, the frost-limited P. arborea. Based on lapse rates, with 65–130 m lower sea levels during 35–18 ka (Lambeck et al., 2014), and today's distribution of P. arborea on TdC and Gough Island (Ryan, 2007), we can estimate that its absence after 28 ka implies minimum winter temperatures at least 3 ^∘C lower than today, which agrees well with our MAAT curve (Fig. 3b).

To evaluate changes in the degree of weathered material, we used a modified Chemical Index of Alteration (CIA) (Fig. 4d). The long-term development can be divided into three phases with initially low but very variable values until 31 ka, a second phase with stable intermediate CIA values until 27 ka, followed by higher and varying values (Fig. 4d). Magnetic susceptibility (MS) shows centennial–millennial oscillations superimposed on an increasing trend from the bottom to the top of the core (Fig. 4m), which is an indicator of in-washed mineral matter from magnetite-rich basaltic rocks of the catchment. The values of total organic carbon (TOC) and biogenic silica (BSi) (Fig. 4l and j) reflect organic and aquatic productivity in and around the lake with highest values in the oldest section. TOC shows a general decline and BSi oscillates with higher values until 28 ka, after which it gradually drops. The fairly high C∕N ratios (Fig. 4k), with a mean value of 17.6, show that organic matter is a mix of terrestrial and aquatic sources. The high and oscillating ratios in the older section followed by a gradual decline implies terrestrial sources dominating until 28 ka, after which time aquatic sources become more important. With respect to bulk stable isotopes (Fig. S1), the high δ¹⁵N values imply a marine influence possibly related to the presence of marine birds (Caut et al., 2012), such as great shearwater and albatrosses, which have a great impact on the ponds today. Rising δ¹³C values at 25.7 ka are consistent with the declining C∕N ratios after 28 ka, i.e., more aquatic material with enriched ¹³C and perhaps in combination with a higher influence from C₄ grasses.

Unlike the pollen record (Fig. 5), the diatom record shows large shifts and the 33 diatom taxa (Fig. S2) have been classified into three environmental forms. Changes in these groups imply shifts in aquatic and environmental conditions in and around the lake. They show a lake with open water early in the record, followed by shifting lake levels between 35 and 33 ka (Fig. 4i), supported by δD values of long- and mid-chain n-alkanes (Fig. 4f). At 31 ka the open-water ratios drop and reach a minimum at 29 ka, in antiphase with the acid water diatom ratios (Fig. 4g), followed by a rise until 26.6 ka. Thereafter acid species dominate as oligotrophic wetland encroached around the lake, while periods of more terrestrial diatoms imply episodes of in-washed diatoms from the surroundings. Around 21.2 ka more open-water conditions prevail again with high ratios during 19–18.6 ka, before the lake dried out (Ljung et al., 2015). The shifts in diatom communities show that 1P went through substantial hydrologic changes, some of which were rapid, induced by changing P∕E ratios, in contrast to the fairly stable vegetation around the lake as seen in the pollen record.

4.2 Linking the Nightingale Island record to South Atlantic hydroclimate

The hydrological sensitivity of a basin like 1P makes it ideal to place local changes into the context of regional hydroclimate shifts. To analyze the variability through time Principal component analysis (PCA) was carried out on a data set with 14 hydroclimate-sensitive proxies resulting in three PCA zones (Fig. 6a). Note that the resolution of the PCA record depends on the proxy with least common sample levels (Table S2), in this case biomarker analyses. Therefore the temporal resolution of the PCA is not as high as some ice core and marine records. Based on the proxy loadings in the PCA plot (Fig. 6b), it can be divided into four different segments with variable hydroclimate and environmental conditions. The importance of presumed temperature proxies on Axis 1 (38.1 % of the variance) is evident where warm pollen ratios, Phylica arborea pollen, BSi, TOC and open-water diatoms show warm humid conditions to the left (negative) in the biplot (Fig. 6b) vs. cooler and drier to the right. The latter is accentuated by Ophioglossum spores, a fern growing at high and cold altitudes on TdC. A correlation analysis between the PC1 values and our MAAT values (Fig. 4b) shows an r² value of 0.56, corroborating that Axis 1 mainly represents temperature. Axis 2 (13.4 % of the variance) is linked to hydrologic indicators being dominated by the δD values of the aquatic n-C₂₁₋₂₃ and terrestrial n-C₂₇₋₃₁ alkanes (Fig. 6b). We interpret higher δD values (positive Axis 2 values) to show stronger influence of more local air masses, with more evaporation and semiarid conditions, also shown by Empetrum rubrum pollen in the upper left quadrant, while the upper right quadrant of the plot shows an acid oligotrophic swampy setting. The segment to the lower right in Fig. 6b displays cold conditions and in-wash of terrestrial diatoms as an effect of higher lake level during episodes of more precipitation. The lower left represents warm and wet conditions, implied by P. arborea pollen and open-water diatoms and, in general, negative Axis 2 values relate to more negative δD values.

To illustrate the relation between the position of the westerlies and the isotopic composition of precipitation at TdC in the simulation, we selected extreme values of high and low δD at TdC (Fig. 7a) and made composite anomalies of the annual mean westerly wind strength at 850 mbar (u850mb, Fig. 8a), precipitation (Fig. 8b), 2 m temperature (t2m, Fig. 8c) and precipitation-weighted δD (Fig. 8d) for high minus low δD at TdC. This shows that the variability in δD in precipitation at TdC is only weakly dependent on local temperature. Instead, shifts in δD at TdC are related to large-scale changes in precipitation and the position of the westerlies. Positive δD anomalies at TdC imply a more southern position of the core of the westerlies with drier and more subtropical conditions at TdC, and negative δD anomalies at TdC denote a more northern position of the core of the westerlies bringing more polar air masses with wetter conditions at TdC. From Figs. 7 and 8, we note that the shifts in TdC precipitation are governed by the precipitation zone on the northern flank of the westerlies shifting with the position of the westerlies themselves. We therefore conclude that our model analysis shows that isotope variability in precipitation at TdC is mainly related to shifts in large-scale circulation. High δD values at TdC imply a more southerly SHW position with stronger winds in its core, while low δD values show a more northerly SHW position with weaker winds (Fig. 7e and h). Our analysis also shows that high (low) δD values are related to less (more) precipitation at TdC but shows little dependency on temperature (Figs. 7f and i, and 8c). The amplitudes of the δD values in our proxies are significantly larger than the modeled amplitudes, implying larger changes in the climate variables during the recorded isotope shifts compared to the year-to-year variability in the model. Furthermore, the modeled relationship between δD and precipitation corresponds well to the PC2 variability in the proxies (Fig. 6b); for example, high PC2 and δD values relate to more Cyperaceae (lake overgrowth) and Empetrum pollen values (arid soils) and more acid diatoms (swampy), while low PC2 values relate to open-water (lake) and terrestrial (flushed-in) diatoms.

Figure 8 Composite maps of changes in wind, precipitation, temperature and δD related to δD variability at TdC, showing annual anomalies based on composites for high and low δD at TdC (see Fig. 7a). (a) Westerly wind speed (850 mbar U wind, positive towards east (m s⁻¹)). The dashed gray lines show the approximate northern and southern boundaries of the westerlies (850 mbar U wind > 5 m s⁻¹) to clarify that high TdC δD is related to a southwards shift in the westerlies. (b) Precipitation (mm month⁻¹). (c) 2 m temperature (t2m (^∘C)). (d) Precipitation-weighted δD (‰). Stippling indicates significant anomalies (p < 0.01) calculated using a two-tailed Student's t test. The green spot shows the position of TdC.

5 Hydroclimate correlations and interpretations

5.1 The large-scale hydroclimate pattern

The three PCA zones, dated to 36.2–31.0, 31.0–26.5 and 26.5–18.6 ka, show a trend and pattern which is recognizable in much of our data set as well as in the EDML (Fig. 4a) and South Atlantic marine record (Fig. 10b). Zone 1 is fairly warm but oscillates between low and high PC2 values, related to more northerly and weaker SHW, and more local air masses with stronger westerlies in a more southern position, respectively. Zone 2 is generally more stable with some minor oscillations with more southerly SHW and corresponds largely to the fairly warm period in Antarctica with the three isotope maxima, AIM4.1, AIM4 and AIM3 (Fig. 4a), and a stable and mild period in the South Atlantic marine realm (Fig. 10b). Zone 3 shows a cooling trend, also visible in the EDML and marine record, with variable SHW. It appears that TdC was continuously influenced by the SHW, as shown by the absence of arid conditions and generally low δD values, verified by humid conditions in southwestern-most Africa throughout most of MIS3 and MIS2 (Chase and Meadows, 2007). Apart from the resemblance between the long-term trends in Antarctic ice core data and marine data at 41^∘ S in the South Atlantic (Barker and Diz, 2014) with our data it is, in spite of our lower resolution, interesting to compare our PC2 and δDn-C_C27−C31 records (Figs. 4a and 10g) with other regional records related to SHW. Taking age uncertainties in a few hundred years into account we note a resemblance with marine Fe fluxes at 42^∘ S (Martínez-García et al., 2014) where low δD values (Fig. 10g) covary with high Fe fluxes (Fig. 10f) due to northerly SHW in a cooler Southern Hemisphere, thus expanding the Patagonian dust source. Similar covariability can be seen in the δ¹⁸O record on fluid inclusions of SE Brazilian speleothems (Millo et al., 2017), where low values (Fig. 10e) imply strengthening of the monsoon shifting the South Atlantic atmospheric system southwards, including SHW. We also note that the Antarctic CO₂ record (Fig. 10c) and the ${\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}$ record (Gottschalk et al., 2015) from the South Atlantic (Fig. 10d), inferring AMOC strength and Southern Ocean ventilation, share similarities with our SHW records, described in the section below.

Figure 9 Parametric plot of the PC1 and PC2 sample values as a function of time shown by the color bar to the right. Red numbers denote each age with gray numbers at the start and end of the plot. Data were interpolated to 50-year time steps to illustrate rate of change; the larger the distance between dots, the more rapid the change. Note that the hydroclimate interpretations from Fig. 6b are shown on the two PC axes.

5.2 A detailed hydroclimate scenario for the central South Atlantic

Due to chronological uncertainties in all records, lower resolution in some records and the complex phase relationships during abrupt interhemispheric climate shifts (Markle et al., 2016), detailed comparison of short-term variations across sites has to be treated with caution. In spite of these short-comings we will present a scenario based on our record and likely correlations.

The start of our record shows warm and wet conditions with northerly SHW, coinciding with the long and warm AIM7 followed by a cooling (Fig. 10j) at the onset of DO7. This is followed by the very dynamic period, 35–33 ka, shown by high sedimentation rates (Fig. 3) and peak variability in terms of both rapidity and amplitude (Fig. 9). Such variability is also seen in marine and ice core records, and in spite of the age uncertainties at 34–35 ka (Fig. 3) we tentatively correlate this period in our record to the end of DO7 and the minimum between AIM6 and AIM7. This corroborates the overlaps and time lags that have been postulated for DO and AIM events (Markle et al., 2016; Pedro et al., 2018; WAIS Divide Project Members, 2015). At 34 ka we note a temperature peak at the onset of AIM6 (Figs. 4 and 11) followed by falling temperatures, δD, Fe flux and CO₂ values and high humidity (Fig. 10g, f, c and h). This change reflects northerly and weaker westerlies, with rising speleothem δ¹⁸O and West Antarctic Ice Sheet (WAIS) d_ln values (Figs. 10e and 11e), denoting the start of DO6 with a warming of the NH (Fig. 10a). This caused a northwards-shifting ITCZ and SHW in line with the theory that the atmospheric circulation system moves towards the warmer hemisphere, responding to the change in the cross-equatorial temperature gradient (McGee et al., 2014). At 33.5 ka we see a southwards SHW shift with rising temperatures and higher CO₂ and lower WAIS d_ln values with dry conditions. We relate this to the onset of AIM5: a warming which is interrupted at 32.8 ka by a northerly SHW shift and wetter conditions (Figs. 11f and 10h) possibly triggered by DO5. This partly continues until 31.7 ka when SHW move south with a minor temperature rise (Fig. 11d) and decreasing humidity, possibly as a response to the post-DO5 cooling (Fig. 10a). The high variability and large amplitude of the changes in Zone 1 (Fig. 9) have facilitated conceivable correlations to other records. Based on these results we can conclude that the variability of PC2 implies a northerly shift of the SHW and a more southerly position during warm periods in Antarctica, also in line with interpretation of Antarctic deuterium excess data (Markle et al., 2016).

Figure 10 Comparisons between other proxy records (a–f) and Nightingale Island proxies for SHW (g), wetness (h) and temperature (i), with mean values as broken lines. (a) δ¹⁸O values from the North Greenland Ice core Project (NGRIP) ice core (Andersen et al., 2006) showing DO and H events. Ice core records are on a common timescale (Veres et al., 2013). (b) Abundance (%) of polar foraminifera at 41^∘ S in the South Atlantic (Barker and Diz, 2014). (c) Composite Antarctic CO₂ record from Siple Dome (Ahn and Brook, 2014) and WAIS (Stenni et al., 2010). (d) (${\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}$) data at 44^∘ S in the South Atlantic (Gottschalk et al., 2015). (e) Speleothem ¹⁸O record on fluid inclusions from SE Brazil (Millo et al., 2017). (f) Fe flux data in the South Atlantic at 42^∘ S (Martínez-García et al., 2014). NI data are then followed. (g) Average δD values for the terrestrial n-C₂₇₋₃₁ alkanes. (h) Abundance (%) of temporary phytoplanktonic diatoms implying relative water depth. (i) MAAT from the GDGT analyses. Note that sample levels are shown by a dot in (g)–(i) and that y axes of (b) and (g) show higher values downwards to facilitate comparisons to other proxies. Note that the two thick gray lines (31 and 26.5 ka) indicate the position of the three PCA zones (Fig. 6).

The Zone 1–Zone 2 boundary at 31 ka (Fig. 6a) is a dynamic transition, shown by many proxies and peak sedimentation rates (Figs. 9 and 3). The 4.5 kyr long and stable Zone 2 (Fig. 6a) is characterized by fairly high but slightly decreasing temperatures and, as in Zone 1, a dominating southerly SHW position. It is possible, taking age uncertainties into account, that H3 at 30.5 ka (Fig. 10a) triggered the southbound SHW, the rising CO₂ and MAAT values, and the reduced humidity between 31 and 30 ka (Figs. 11f, 10c, 10i and 10h). The following long and warm AIM4 may have stabilized conditions in the South Atlantic in spite of the DO4 event at 28.8 ka. This stability is also seen in marine records (Fig. 10b), and the rather stable southern position of the SHW agrees with the fairly high CO₂ values between 30 and 27.2 ka and with falling and rather low Fe fluxes (Fig. 10f). We also note higher lake evaporation from δD values of the aquatic n-C₂₃ (Fig. 4f) in concert with rising summer insolation (Fig. 11a). Around 27.5 ka we see a brief response in some of the proxies to the short DO3 event (Fig. 10a), such as the MAAT and PC2 records (Fig. 11d and f), which is also noticeable in, for example, the marine and Brazilian monsoon records (Fig. 10b and e).

The start of Zone 3 at 26.5 ka (727 cm) constitutes the most drastic change in our record (Figs. 6a and 9) but timing varies between proxies (Fig. 4). MAAT, TOC and C∕N ratios start to decrease already at 28–27.5 ka, coinciding with DO3, while the biologic proxies (Figs. 4C and g–j, S2) respond slightly later possibly because they do not react until certain hydroclimate thresholds for the vegetation and algae flora are reached. The Zone 2–Zone 3 transition is roughly simultaneous with the onset of LGM in Antarctica (Fig. 10c), when 1P switched from a lake to a wetland, coinciding with increased abundance of polar foraminifera at 41^∘ S (Fig. 10b). This may be an effect of the SF moving north of TdC, a meridional shift comparable to what has been shown from the eastern Pacific (Kaiser et al., 2005). The fairly stable PC1 values show cool and less-humid LGM conditions, while the variable PC2 values imply shifts in the position of SHW (Fig. 11f). There is also some correspondence between our δD (n-C₂₇₋₃₁) maxima after 27 ka and Fe flux minima from the South Atlantic (Fig. 10g–f), both indicating southerly shifts of SHW. During this period our data also show generally higher mean δD (n-C₂₇₋₃₁) values than in Zone 1, implying a more southern position of SHW during the Antarctic LGM, as seen in some modeling results (e.g., Sime et al., 2016). This is also compatible with the fact that the LGM temperature lowering in the Northern Hemisphere (Johnsen et al., 1995) was much larger than in the south (Stenni et al., 2010), shifting the atmospheric system to the south due to changes in the cross-equatorial gradient (McGee et al., 2014), also implied by the speleothem δ¹⁸O data (Fig. 10e) showing increased precipitation (Millo et al., 2017).

Figure 11 Comparison between our PC1 and PC2 records and other relevant data. (a, b) Mean daily summer and winter insolation at 37^∘ S (Laskar et al., 2004). (c) EDML δ¹⁸O record (EPICA Community Members et al., 2006) with Antarctic Isotope Maxima (AIM). (d) PC1 scores implying temperature shifts at NI. (e) WAIS d_ln values from west Antarctica (Markle et al., 2016). (f) PC2 scores indicate impact of SHW at NI. Note that sample levels, i.e., time resolution, for the PC records are shown as dots. Note the two thick gray lines (31 and 26.5 ka) indicating the position of the three PCA zones (Fig. 6).

After 26.5 ka we note phases of less-humid swampy oligotrophic conditions on NI at 26, 24.5–23, 22 and 20.5–19 ka (Fig. 10h) interrupted by periods of more or less open water, possibly driven by shifts of SHW. The former often show enriched δD values (Fig. 10g), while the latter were characterized by higher precipitation and more depleted δD values. Regarding the response of CO₂ to these SHW shifts, we note a fairly good agreement between low/falling CO₂ values and a northerly SHW position, and vice versa. For example, the CO₂ minimum at 24.5–25 ka (Fig. 10c) matches with an extreme northern SHW position (Figs. 10g and 11f), and the CO₂ peak at 23.3 ka agrees with the end of a long phase of southwards-moving SHW. The latter might have been triggered by the onset of H2 at 24.1 ka (Fig. 10a) followed by the inception of AIM2 (Fig. 11c).

The absence of P. arborea (Figs. 4c and 5) and our temperature proxies (Fig. 4b) implies that minimum winter temperatures at our site were occasionally below zero, especially after 26 ka; periods of frost also explain increased mechanical weathering (Fig. 4d). Between 23 and 19 ka the Antarctic winter sea ice reached 47^∘ S in the South Atlantic (Gersonde et al., 2005), only some 1000 km south of TdC. Our 1P record shows a declining temperature trend during the end of this period (Fig. 10i), in contrast to rising temperatures in Antarctica and the South Atlantic (Figs. 11c and 10b). This regional temperature anomaly may be explained by the declining summer insolation at the latitude of Tristan da Cunha (Fig. 11a) and may also, at the end of the LGM, be related to the break up of Antarctic ice shelves as sea levels rose, causing cooler conditions further north. In fact, temperature minima after 19 ka are seen in both our record and in marine data (Fig. 10i and b), as well as a δD minimum (Fig. 10g).

5.3 A climate synthesis

In general, our data imply two main climate modes for the study period, separated by a transition period of 31–26.5 kyr, Zone 2. This is displayed in Fig. 9 with pre-LGM (Zone 1) clearly separated from the LGM period (Zone 3) on Axis 1 but also with higher variability in the pre-LGM period. This variability is possibly related to an active bipolar seesaw mechanism during Zone1 or MIS3 even at the fairly low latitudes of TdC, triggering N–S shifts of SHW and related hydroclimate conditions. Any CO₂ effects from the rapid SHW shifts in Zone 1 are not discernible but the dominating, more northern SHW position may have resulted in the general CO₂ decline (Fig. 10c). With the onset of Zone 2 there may be a stronger link between CO₂ and SHW. In view of carbon cycle time lags, the mainly southerly-positioned and more intense SHW at 31–27.5 ka (Fig. 11f) may have resulted in the rising and higher CO₂ concentrations at 30.5–27.2 ka (Fig. 10c), with more upwelling, CO₂ outgassing and less sea ice. The LGM mode is characterized by falling and low temperatures and a lack of clear effects of the bipolar seesaw mechanism possibly due to the much stronger cooling in the north as the cross-equatorial gradient changed. The variability is mainly related to proxies associated with SHW changes, as summarized by PC2, with a similarly high-frequency variability in WAIS d_ln and Fe fluxes (Figs. 11e, 10f), with resulting CO₂ variability. However, a key difference between our SHW proxies (PC2) and the WAIS d_ln record is that the latter represents SHW variability superimposed on large-scale temperature trends while our PC2 record reflects the SHW signal without temperature impact.

Thus, the largest change in our record occurs after 27.5 ka when the effects of the strong post-DO3 cooling of the Northern Hemisphere start dominating the hydroclimate of the South Atlantic with highly variable SHW after 25 ka; possibly a prerequisite for the oscillating CO₂ levels after the CO₂ minimum at 25 ka (Fig. 10c).

6 Conclusions

In conclusion, we think our Nightingale Island data demonstrate the potential for remote-island proxy records to register large-scale atmospheric shifts in an oceanic setting, especially if the island location, in relation to marine and atmospheric fronts, is chosen well. In addition, if the right types of proxies are chosen, multiproxy lake records have a particularly large potential since they disclose both terrestrial and aquatic responses to shifting atmospheric conditions, as shown in our 1P record. By combining these responses they can be translated into relative changes in hydroclimate conditions.

Our 1P (1st Pond) data show that the glacial hydroclimate of South Atlantic mid-latitudes experienced varying degrees of humidity but with a more or less continuous impact of the SHW. Temperature conditions were in general warm but oscillating during MIS3, with shifting strength and positions of the westerlies. Weaker and northwards-moving SHW at the onset of NH interstadials with stronger and southerly westerlies during NH stadials partly reflect the complex processes behind phase relationships between Greenland and Antarctic ice core climate records (Pedro et al., 2018). These shifts, possibly triggered by changes in the cross-equatorial gradient, are to some extent manifested by rising (falling) CO₂ levels when SHW were stronger (weaker) and located more towards the south (north), in line with Holocene records (Saunders et al., 2018). The largest variability in our record is seen during the fairly warm and humid period of 36.5–31 ka with frequent and abrupt shifts, followed by a fairly stable period, 31–27 ka, with slowly declining temperatures and dominating southerly SHW. The largest overall change occurs after 27 ka, exhibited by a distinct cooling trend. This early mid-latitude cooling is in phase with LGM in Antarctica, consistent with some modeling results (Fogwill et al., 2015). We think this represents a mode shift in hydroclimate: from the highly variable MIS3 conditions through the more steady conditions during 31–27 ka (Fig. 11d and f) into LGM with its cool and less-humid climate, perhaps as a result of the SF moving north of TdC. The variable position of SHW (Fig. 11f), with particularly high δD values at 26, 23.1, 21 and 19.1 ka (Fig. 10g), is noteworthy, inferring fairly sudden and distinct southerly shifts of the westerlies. The end of our record shows that fairly cool conditions persisted in these SH mid-latitudes until at least 19 ka. This might have been a combined effect of declining summer insolation and northwards-shifting westerlies (Fig. 11a and f), conveying cold air masses, sea ice and ice bergs far north of the first peak of iceberg-rafted debris by Weber et al. (2014) from the breaking up of Antarctic ice shelves starting at 20 ka, and named MWP-19KA.