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 3Age 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.

Download

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 4Antarctic 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).

Download

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 5Pollen 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.

Download

Figure 6Principal 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).

Download

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 7Zonal 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.

Download

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 8Composite 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.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 9Parametric 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.

Download

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 10Comparisons 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).

Download

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 11Comparison 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).

Download

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).