Shut down of the South American summer monsoon during the penultimate glacial.

We analysed changes in mean annual air temperature (MAAT), vegetation and biomass burning on a long and continuous lake-peat sediment record from the Colônia basin, southeastern Brazil, examining the responses of a wet tropical rainforest over the last 180 ka. Stronger southern atmospheric circulation up to the latitude of Colônia was found for the penultimate glacial with lower temperatures than during the last glacial, while strengthening of the South American summer monsoon (SASM) circulation started during the last interglacial and progressively enhanced a longer wet summer season from 95 ka until the present. Past MAAT variations and fire history were possibly modulated by eccentricity, although with signatures which differ in average and in amplitude between the last 180 ka. Vegetation responses were driven by the interplay between the SASM and southern circulation linked to Antarctic ice volume, inferred by the presence of a cool mixed evergreen forest from 180 to 45 ka progressively replaced by a rainforest. We report cooler temperatures during the marine isotope stage 3 (MIS 3: 57-29 ka) than during the Last Glacial Maximum (LGM: 23-19 ka). Our findings show that tropical forest dynamics display different patterns than mid-latitude during the last 180 ka.

Introduction

Tropical rainforests (e.g. Amazon, Atlantic forest) cover 49% of the South American territory[1] but are of global importance as they account for nearly 50% of total carbon storage[2] and the water cycle[3,4], and host between half and two-thirds of the world’s species[5]. The Atlantic forest domain offers a unique opportunity to test palaeoclimatic and palaeoenvironmental dynamics in tropical rainforests, because its floristic composition has been determined by fluctuating climatic conditions in an altitudinal gradient[6] (from 0 to 2,900 m a.s.l.). The forests are located on a long latitudinal gradient along the Atlantic Ocean coast (0° to 23° S) where the intertropical convergence zone (ITCZ) and the SASM control the precipitation regimes (Fig. 1). Recent studies showed that long-term spatial variations of the hydroclimate led to different SASM responses in the Neotropics[7,8]. In Southeastern Brazil, the SASM appears to be paced by a precession signal[9,10]. However, the question whether the tropical rainforests in South America expanded, retreated or changed their floristic composition during glacial periods remains mostly unanswered due to data scarcity and model discrepancy[11-13]. Understanding ecosystem responses to climatic transitions such as those observed during glacial cycles is important when considering the long-term implications of the current climatic crisis (e.g. rapid increase in temperature, sea level and atmospheric CO2), because the future of tropical ecosystems indeed depends on how fast they adapt to the different drivers of change over time.

Figure 1

Regional setting and study site. (a) Map of South America, showing the distribution of precipitation, plotted as the long-term mean annual precipitation from 1981 to 2010 using the data from the Global Precipitation Climatology Project (GPCP)[85]. Highlighting the South American summer monsoon (SASM; black arrows), and the position of the Intertropical and South Atlantic convergence zones (ITCZ and SACZ; blue dashed lines) during the austral summer in South America (DJF). Blue arrows show the northernmost position of the year-round southern fronts. Mean monthly precipitation at Colônia[29] (black bars) and mean monthly temperature in São Paulo[30] (red line) are shown at the bottom of the panel. The numbers refer to the regional datasets to be compared with Colônia (yellow star); (1) Botuverá[9] (2) Santana[42] (3) Sea surface temperature record GL-1090[39]. (b) Map showing the distribution of the Atlantic forest in Brazil. (c) Geology of the Colônia structure showing the location of the sediment cores CO3 and CO14 (adapted from Riccomini et al.[24]).

So far, only four lowland rainforest pollen records covering more than one glacial-interglacial cycle exist. In northern Australia and West Africa, the general forest retreat corresponded to a worldwide decrease in mean annual precipitation during glacial periods[14,15]. In northern Amazonia (Lake Pata), the rainforest continued throughout the last glacial-interglacial with no significant difference in floristic composition between marine isotopic stage 6 (MIS 6: 191 to 135 ka) and MIS 5 (135 to 71 ka)[16]. The Southern Hemisphere records (Lynch crater[14] and Colônia[17]) show that mixed evergreen rainforests expanded and were more continuous during the penultimate glacial period (MIS 6) than during the last glacial period (MIS 4, 3, 2: 71 to 14 ka). Contrarily, northern mid-latitude records show forest expansion phases occurring only during interglacials[18]. It has been suggested that the Antarctic ice volume and atmospheric pCO2 modulate the expression of glacial periods in the Southern Hemisphere, particularly the dynamics of the subtropical front and ITCZ[19,20]. This highly dynamic system supplies moisture for tropical precipitation in a round-trip circuit driven by obliquity cycles[21,22]. However, if obliquity drives Southern Hemisphere climate and Antarctic sea ice expansion[23], the tropical precipitation cycle appears to be mainly driven by precession[9]. One can therefore wonder how the southern rainforest responded to differences in the expression of glacial severity between palaeoclimatic forcings such as temperature and hydroclimate.

Here, we explore the responses of a wet tropical forest in the Atlantic forest domain to changes in MAAT and precipitation during the last 180 ka. We compare the observed regional responses with Southern Hemisphere climate records (Antarctic ice volume and pCO2) to characterise their respective behaviours.

Study area

The Colônia basin (23°52′03″S and 46°42′27″W, ca. 700 m a.s.l.) is a depression surrounded by an annular rim with hills reaching 125 m (ca. 900 m a.s.l) (Fig. 1). The crater-like structure was formed in crystalline basement rocks, mostly granitic gneiss, some schist and quartzite of Neoproterozoic age (600–700 Ma), presumably by a meteorite impact[24,25]. The inner structure is filled by alluvial sediments capped by a swamp partially drained by the Vargem Grande stream, which flows eastward through the eastern rim[24] (Fig. 1).

The regional climate is strongly influenced by northward shifts of polar air masses that result in permanent drizzle and moisture in winter[26,27] (JJA), while during austral summer (DJF) it is affected by the SASM and the position of the South Atlantic Convergence Zone (SACZ)[28] (Fig. 1). Mean annual precipitation is 1,600 mm, including two months with less rainfall in austral winter (a short dry season in July and August), although natural fires do not occur during this period due to high moisture availability in the area. The mean annual temperature is 20 °C and the average winter temperature is ~17 °C (Fig. 1)[29,30]. The long-term (from 1981 to 2010) mean monthly temperature anomalies between the city of São Paulo[30] (where Colônia basin is located) and sea surface temperature (SST) near the coast[31], reflect less negative differences (−0.2 ± 0.6) during the peak of the influence of the SASM (DJF). The anomalies become more negative (−6 ± 0.5) during winter (JJA) when the precipitation is lower in the region (weaker SASM).

The rim surrounding the Colônia basin is covered by Atlantic rainforest mainly represented by the families Myrtaceae, Rubiaceae, Bromeliaceae, Podocarpaceae, Arecaceae, along with other evergreen forest species. The swamp located in the centre of the basin is dominated by Poaceae, Xyridaceae (Xyris), Lentibulariaceae, Sphagnaceae (Sphagnum) and Asteraceae[32].

Results

In 2014, a 14-meter sediment record (CO14) was recovered in the swampy area of Colônia using a D-section Russian corer[25]. A sharp sedimentological change at 850 cm revealed a transition from a lacustrine to a peat environment (Supplementary Fig. S1). Distinct sedimentary deposit regimes were characterised and differentiated using X-Ray fluorescence measurements (XRF) (Supplementary Fig. S9). Biological proxies including pollen, branched glycerol dialkyl glycerol tetraethers (brGDGTs) and microcharcoal content were used to reconstruct the vegetation, MAAT and biomass burning. The Colônia composite pollen record contains data from two parallel sections (CO3 and CO14). CO3[17] was used to compare new proxies of CO14 from mid-MIS 5 until present (see Methods and Supplementary Figs. S5 and S6). Pollen interpretations are based on the arboreal and non-arboreal pollen curve, cool and moist taxa such as Araucaria and Podocarpus, the most relevant evergreen forest taxa indicators (Anacardiaceae, Ilex, Myrsine, Myrtaceae) and semi-deciduous taxa indicators (Alchornea, Celtis, Euphorbiaceae, Malpighiaceae and Melastomataceae)[33,34]. MAAT reconstruction is based on the brGDGTs temperature calibration from Naafs et al.[35], with temperature values ranging from 13.5 to 26 °C (see Methods section and Supplementary Fig. S7 for a discussion on the GDGTs sources and calibration at Colônia). Biomass burning interpretations are based on microcharcoal concentration data expressed as the total number of microcharcoal particles per gram of dry sediment (nb.g−1) (see Methods section). The chronological framework was obtained using a Bayesian approach[36] (See Methods section and Supplementary Fig. S1). It relies on 22 radiocarbon dates (Supplementary Table S1), the Laschamp[37] geomagnetic excursion and one optically stimulated luminescence (OSL) date as a validation point for the chronology (Supplementary Table S2). Additionally, the obtained MAAT was fitted to the benthic δ18O LR04 stack[38] using three changes in temperature as tie points (See Methods section, Supplementary Fig. S2, Table S3). The obtained chronology spans the last 180 ka, with a very well defined age-depth model at the top of the core (approximate mean uncertainty of ±1.1 ka) and a loose model below the 14C limit (approximate mean uncertainty of ±12 ka) (Supplementary Fig. S1).

Changes from glacial to interglacial are reflected in the temperature brGDGTs-derived curve (MAAT), biomass burning and the first principal component (PC1) calculated using the XRF core scanning dataset (PC1 accounts for 79% of total variance) (Fig. 2 and Supplementary Fig. S9). Land-sea temperature gradient between the anomalies from the MAAT and the SST GL-1090 record[39] (See Methods section and Supplementary Fig. S8) reflects sustained changes between MIS 5, MIS 2-3-4 and MIS 6. Negative numbers imply relative cooling on land and thus less moisture transport in the region (Fig. 3).

Figure 2

Summary of palaeoecological proxies from Colônia. (a) Orbital parameters for latitude 23°52′03″S using La2004[86]. (b) PC1 of the intensity of the XRF elements from CO14. (c) Reconstructed MAAT (CO14). (d) Biomass burning from CO14. (e) Arboreal (AP) pollen percentages from CO3 and CO14. (f) Podocarpus pollen percentages from CO3 and CO14; *reflects a peak of ca. 60% (see Methods section). (g) Araucaria pollen percentages from CO3 and CO14. The marine isotopic stage boundaries ages are based on the LR04 δ18O stack[38,87]. The boundary between MIS 5 and 6 is based on Henderson and Slowey[88].

Figure 3

Regional comparison during the last 180 ka. (a) Speleothem records (Botuverá (black) and Santana (gray))[9,42] plotted with insolation (W/m2) at latitude 23°52′03″S using La2004[86] (red dashed line). (b) Arboreal pollen percentages (AP) from Colônia (CO3 and CO14). (c) Percentage curve of the sum of cool moist conifers Araucaria and Podocarpus; the black star locates the last continuous occurrence of Araucaria, and asterisk represents ca. 60% of Podocarpus (see Methods section). (d) Land-sea temperature gradient. (e) pCO2 (Lüthi et al.[89]). (f) Total ice volume represented in m s.e. (meters sea-level equivalent) for the Antarctic ice sheet[45]. The marine isotopic stage boundaries ages are based on the LR04 δ18O stack[38,87]. The boundary between MIS 5 and 6 is based on Henderson and Slowey[88]. Pink bands highlight the maximum insolation peaks of influence of the SASM system during the last 95 ka. The blue band also highlights the maximum insolation peak but emphasises the differences in Colônia’s proxies compared with the regional datasets.

The penultimate glacial (MIS 6) was characterised by the coolest temperatures with an average MAAT of 15.6 ± 1 °C and full expansion of the cool mixed evergreen forest (40–70%). The forest was dominated by Myrsine, Myrtaceae, Araucaria and Podocarpus, and semi-deciduous taxa including Melastomataceae, Malpighiaceae and Celtis (Fig. 2 and Supplementary Fig. S5). Open vegetation (30–60%) was dominated mainly by Poaceae and Asteraceae. During MIS 6, biomass burning was on average 11.5 ± 8 × 104 nb.g−1. In comparison, the last glacial period (MIS 4, 3, 2) was warmer with an average temperature of 19.5 ± 2.2 °C and was characterised by gradual temperature changes. During that period, MAAT progressively decreased from 24 °C to a minimum temperature of 15.6 °C at 43 ka, and then increased again up to 22 °C from 41 to 21 ka, meaning temperatures were lower during the MIS 3 than during the LGM[40]. Two fire episodes occurred during the last glacial period with contrasted biomass burning from ca. 71 to 50 ka (16 ± 17.7 × 104 nb.g−1 on average) and from ca. 35 to 23 ka (66.5 ± 66.3 × 104 nb.g−1 on average). Additionally, the vegetation was characterised by fluctuations from mixed evergreen forest to open areas dominated by grasslands[17]. During MIS 4 (71–57 ka) the landscape was mainly characterised by grasslands (9–66%) although a mixed evergreen forest persisted. MIS 3 (57–29 ka) was characterised by early expansion of grasslands and an increase in mixed evergreen forest between 50 and 40 ka, followed by a decline in forest until 30 ka. During MIS 2 (29–14 ka), the mixed evergreen forest reached full expansion between 30 and 25 ka, and grasslands were dominant during the late glacial. Variations in detrital inputs detected by the PC1 from the XRF core scanning (e.g. Si,Ti, inc/coh) suggest increased erosion rates during the dry periods reported by speleothem records and forest contractions over the last glacial period (see description of the sediment in Supplementary Material).

Discussion

Palaeoclimate and vegetation responses to glacial-interglacial cycles

During the penultimate glacial, cooler temperatures in Colônia (on average 4 °C lower than during the last glacial) favoured the presence of permanent cool mixed evergreen forest with conifers Araucaria and Podocarpus associated with strong southern circulation towards the tropics. This has been assumed since recent estimates based on modern pollen climate calibration in the Atlantic forest domain showed that the continuous strong northward shifts of southern air masses were linked to expansion of the cool mixed forest[27,34,41]. In addition, a contrasting land-sea anomalies gradient (−7.4 ± 1.2) could explain a reduced SASM intensity to the study site, due to the effect of cooler land temperatures compared to the SST during that period (Fig. 3 and Supplementary Fig. S8). Such dry effect can be observed in the relatively frequent fires, the high frequency of open vegetation and semi-deciduous forest taxa (Supplementary Fig. S5). On the other hand, the presence of the mixed evergreen forest points out that there was still high moisture input (Fig. 2, Supplementary Fig. S5). Today, natural fires do not occur in the area, not even during the dry winter season (July-August), because of the permanent moisture rates (Fig. 1). We therefore suggest that during the MIS 6, SASM circulation was weaker than today, or nearly absent, implying a summer dry season. This hypothesis relies on observations from two nearby speleothem records showing that the strength of the SASM over the last 112 ka is related to insolation in austral summer at precession scale[9] (Fig. 3).

Temperatures at Colônia were warmer during the MIS 5 (130–71 ka) (increase of 7 °C) than during the previous glacial period (MIS 6) at 23.1 ± 1.8 °C on average. Evidence for the highest continuous expansion of the mixed evergreen forest is the presence of cool moist taxa especially Araucaria between 125 and 95 ka (Fig. 2). In this context, although forest expansion would be expected in the high summer insolation phase[17] (more precipitation as inferred from the speleothem records), the abundance of the cool taxa Araucaria and Podocarpus suggests that other factors prevailed during the warm last interglacial period (Figs. 2 and 3). This pattern was probably linked to enhanced extra tropical circulation and abundant winter precipitation, as revealed by the Santana speleothem record from 130 to 112 ka[42]. On the other hand, the dominance of austral summer insolation on precipitation from 95 ka as also revealed by the speleothem records[9], significantly reduced biomass burning (from 95 to 80 ka: 0.9 ± 0.8 × 104 nbg−1 on average) and increased tree pollen frequencies. From early to mid- MIS 5, land-sea temperature differences were gradually reduced (−0.8 ± 2.4), which is likely reflecting a stronger SASM influence during summer in the area[43] (Fig. 3).

The last glacial (MIS 2, 3 and 4) was associated with different global changes in ice volume, temperature and sea level[44,45]. An open landscape associated with drier climate conditions and two phases of erosion (high PC1) characterised the late glacial during low precession phases (with weak SASM activity[9]) and no extra tropical circulation reached Colônia’s latitude. Additionally, the retreat of the cool and moist taxa indicators points to reduced moisture availability in winter. Similarly, the temperature gradient between land and sea increased (−4 ± 2), likely reflecting a decrease on moisture input due to a weak SASM activity. The first phase of biomass burning (71 to 50 ka) suggests that fires were caused by a change in seasonality with an increase in precipitation during austral winter and a decrease during the austral summer. Additionally, the presence of a cool mixed evergreen forest in the region with Araucaria and Podocarpus[34,41,46] points to the influence of southern air masses in southeastern Brazil. The strengthening of the SASM and high summer precipitation during a high precession phase[9] between 50 and 40 ka reduced biomass burning and increased the mixed evergreen forest cover (Fig. 2, Supplementary Fig. S5), also corresponding to bigger differences between land and sea temperatures. The last continuous occurrence of Araucaria observed at 45 ka shows the declining influence of the extra tropical circulation and the increasing dominance of the SASM system in southeastern Brazil. The second episode of biomass burning (35 to 23 ka) in this glacial period is the biggest reported in the record. It started with a low precession phase and a weaker SASM (drier conditions), which does not appear to be linked to grassland or forest expansion (30–25 ka) (Fig. 2). Temperatures at Colônia gradually increased from 15.7 to 22 °C between mid-MIS 3 and mid-MIS 2. Neither the Colônia nor the SST records (GL-1090)[39] contain evidence for cooling during the Last Glacial Maximum, but earlier during the MIS 3.

The relatively cooler Early Holocene, with a slight decrease of 1 °C in MAAT, was followed by a rise in temperature from 3.6 ka on. Renewed biomass burning from 2 ka on, including the highest peak of the record (178 × 104 nb.g−1), coincided with forest regression between 2.3 and 1.6 ka.

Effect of the Southern Hemisphere on the Tropics

The observed environmental and climatic changes between glacial-interglacial periods in Colônia could be orbitally driven at eccentricity scale although a longer record is needed to further test our hypothesis (Fig. 2, Supplementary Fig. S10). Additionally, due to the latitudinal location[47], rainforest expansion and contraction phases were a unique response as a result of a combination of winter seasonality driven by obliquity along with the SASM driven by precession (Figs. 2 and 3, Supplementary Fig. S10). However, additional mechanisms other than insolation cycles must be inferred to explain changes in fire activity or the continuous presence of a cool mixed forest from the penultimate glacial period (MIS 6) until the mid-MIS 5.

The timing and climatic variations observed in the CO14 record could be linked to the differences between temperature and precipitation observed at different latitudes as suggested by models[48]. We propose that the influence of variations in the Antarctic ice sheet on the tropical hydrological cycle affected the floristic composition of the rainforest from the penultimate glacial to mid-MIS 5 (180 to 95 ka). Based on studies of the last glacial in the Southern Hemisphere[44,49,50], we propose that the combination of low PCO2, extended Antarctic ice volume, and an equatorward shift of the westerlies prevented the tropical circulation from influencing the Colônia site from MIS 6 to mid-MIS 5[49]. Those conditions would have allowed the development of a mixed cool evergreen and semi-deciduous forest with no modern analogue[34], along with fire activity during a drier summer season (Fig. 3, Supplementary Fig. S5).

At the beginning of MIS 5, a poleward shift of the westerly belt, an increase on pCO2 and upwelling in the south Atlantic enhanced tropical circulation during austral summer[44,50], although not enough to remove winter precipitation inputs[27,51-53] in the Colônia area (inferred by the presence of Araucaria and Podocarpus). Thus, the Antarctic ice sheet configuration, which controls the position of the subtropical front[19], would explain the shift from dominance of extra tropical circulation[27] during 180 to 95 ka to dominance of the SASM system on the climate patterns from 95 ka onwards (Fig. 3). Throughout MIS 5, during Antarctic ice sheet contractions, the intensity and amplitude of winter rainfall also progressively decreased and a transition phase (80 to 60 ka) characterised by broad changes in forest composition and reactivation of biomass burning, is explained by oscillations in the interplay of the southern air masses and the SASM system.

During the last glacial, the retreat of the cool mixed evergreen forest coincided with a gradual increase in temperatures attributed to the increasing influence of the tropical circulation paced by the precession signal (Fig. 3), while the progressive rise in temperature over the 41 to 21 ka interval was antiphased with low CO2 values and high Antarctica ice volume. Our research highlights the teleconnections between low and high latitude regions and the growing influence of the SASM circulation on climate and vegetation in the southern tropics. The changes in seasonal rainfall distribution from winter centred with continuous fire activity and low temperatures during the penultimate glacial, to summer centred rainfall distribution during the last glacial with major fluctuations in fire activity and temperature, shaped the composition and distribution of the current tropical forest. Further analyses of a deeper core from Colônia are expected to test the long-term eccentricity pacing on regional environmental and climatic changes, the SASM behaviour during previous glacial stages, and characterise how Antarctic ice volume influenced tropical climate and forest diversity.

Methods

Magnetostratigraphy

The half core sections were sampled using u-channels (2 × 2 × 100 cm) for magnetic measurements. Natural, anhysteretic and isothermal remanent magnetisation (NRM, ARM and IRM) were measured at 2-cm depth resolution using a three-axis DC SQUIDS cryogenic magnetometer (2 G 760 SRM) installed in a magnetically shielded room at CEREGE. Alternating field (AF) demagnetisation was applied stepwise at 5, 10, 15, 20, 30, 40, 50, 60 and 80 mT. ARM was imparted in an 80 mT AF and a 0.1 mT direct field (DF) and demagnetised at the same steps as the NRM. IRM was imparted by passing u-channels though a 0.6 T Halbach cylinder to obtain a saturation IRM for magnetite grains (SIRM), and then demagnetised stepwise from 5 to 120 mT. NRM, ARM and IRM display strong stability during AF treatment likely associated with the presence of high coercivity minerals and the small magnetic grain size fraction (single- to pseudo single-domain Ti-magnetite). The demagnetising curves reveal stable intensity, with medium destructive field (MDF) of the NRM barely determined, and stable direction components. Data from the ends of the core (from two to four measurement points) were cleaned to reduce the edge effect associated with the pickup coils response function. An average inclination shallowing of 24° compared to the geocentric axial dipole (GAD) configuration (−12° instead of −36°) likely results from compaction of the organic rich sediments. NRM intensity measured at the 30 mT AF demagnetisation step was divided by ARM intensity to obtain a relative palaeointensity index (RPI). The largest inclination deviation, which reached −87° at a depth of 1.72 m, corresponds to low RPI values and was attributed to the Laschamp[37] excursion, which also fits independent radiocarbon dating. Two other significant inclination deviations associated with low RPI are observed at 3.5 and 6.7 m, respectively.

Radiocarbon dating

A total of 22 samples of bulk peat sediment were analysed at Artemis LMC14 (Gif-sur-Yvette, France) (Supplementary Table S1). The samples were examined under binoculars to exclude possible macro contaminants and were pre-treated chemically using HCl (0.5 N) and NaOH (0.1 N) at 80 °C. Reported 14C activity is expressed as pMC (percent modern carbon), normalised to δ13C of −25‰. The reported ages were calibrated using the IntCal13 curve[54]. The radiocarbon age of sample SacA41593 was calibrated using a post-bomb curve for the Southern Hemisphere[55].

Calibrations and age-depth modelling were performed using R[56] and the RBacon package[36].

Luminescence dating

Luminescence dating was performed at the Luminescence and Gamma Spectrometry Laboratory (LEGaL – Instituto de Geociências, Universidade de São Paulo, Brazil) (Supplementary Table S2). Using opaque syringes, two samples were collected from sandy intervals located at depths of 1 m and 11 m. The surface of the sample that had been exposed to light (ca. 1 cm) was removed under subdued amber light conditions and a total of 20.3 cm3 sediment volume per sample was used for luminescence dating. Grains of quartz and K-feldspar in the 180–250 μm grain size fraction were isolated using standard procedures[57]. Multigrain aliquots of quartz and K-feldspar were measured in a Risø TL/OSL DA-20 reader equipped with blue and infrared (IR) LEDs for light stimulation and a 90Sr/90Y beta radiation source (0.084 Gy/s dose rate). First, equivalent doses were estimated using the optically stimulated luminescence (OSL) signal of quartz[58] detected in the ultraviolet band. The natural OSL signal of the quartz aliquots taken from sample L0247 (1 m depth) was below saturation, but sample L0250 showed saturated natural OSL signals. Consequently, the sediment burial ages of L0250 were determined using the 310 °C isothermal thermoluminescence (ITL) signal of quartz[59] and post-infrared stimulated luminescence (pIRIR) at 290 °C (Buylaert et al.[60]) of K-feldspar (detection in the blue light window). Sample L0250 also showed a saturated natural pIRIR signal, with minimum equivalent dose (2D0) of ~408 Gy. The pIRIR signal consequently only gave a minimum age for sample L0250, and its burial age was determined using the ITL signal from quartz, which was below saturation. We calculated a minimum age based on the IRSL at 50 °C signal (IR50) measured in the first step of the pIRIR at 290 °C dating protocol to constrain the quartz ITL age. The IR50 age is given without correction for athermal fading. The performance of all luminescence dating protocols used in this study was checked through dose recovery tests applied to quartz and K-feldspar aliquots bleached under a solar simulator lamp. The dose recovery test for the quartz OSL signal was performed using quartz aliquots from an adjacent core because sample L0247 did not have sufficient quartz in the target grain-size fraction. The quartz from sample L0247 was kept to estimate equivalent doses. The calculated-to-given dose ratio for the quartz OSL signal was obtained using a preheating temperature of 200 °C and the given dose of 150 Gy was 0.99 ± 0.07 (four aliquots). Given-to-calculated dose ratios for the pIRIR at 290 °C and ITL signals were respectively 1.03 and 1.23 for a given dose of 300 Gy. Equivalent doses of each sample were calculated according to the Central Age Model[61]. Equivalent doses for the sample were calculated with aliquots following acceptance criteria to check corrections for sensitivity changes (recycling ratio within 1.0 ± 0.1) and thermal transfer (recuperation less than 5%). We underline the fact that equivalent doses were calculated using four to six aliquots because of the reduced volume of the available sample.

To estimate the radiation dose rate, natural radionuclide (238U, 232Th and 40K) concentrations were measured by high-resolution gamma ray spectrometry in a high purity germanium detector (HPGe, relative efficiency of 55% and energy resolution of 2.1. KeV) mounted on an ultra-low background shield (Canberra Instruments). Samples (L0250:11 g; L0247:17.4 g) were measured after storage for at least 21 days in sealed plastic containers for radon equilibrium. U, Th and K concentrations were converted into dose rates[62]. The contribution of cosmic radiation to the dose rate was calculated using conversion factors based on the latitude, longitude, altitude, and depth of the samples[63].

The ITL L0250 date was not used in the model and chronological framework due to its large uncertainty (Supplementary Table S2).

Colônia chronology

To check the overall consistency of the Colônia chronological framework, we used three independent posterior validation points:

The OSL age of 31.4 ± 4.3 ka agrees with radiocarbon dating and fits the Bayesian envelope;

we fine-tuned the RPI curve to two independent palaeointensity reference curves[64,65] assuming that the three inclination deviation/low RPI levels record the Laschamp[37] (41 ka) (supported by the radiocarbon dated section in CO14), Norwegian-Greenland Sea[66] (64 ka) and post-Blake geomagnetic excursions[67] (100 ka). Despite legitimate doubts concerning the RPI signal due to complex magnetic mineralogy which prevented its use for straightforward reliable RPI wiggle-matching, the visual relationship within the peat section is acceptable (0–9 m; Supplementary Fig. S3) and support the age model we obtained, i.e., 100 ka at 6.7 m according to the RPI matching versus 109 ± 11 ka derived from the independent chronology;

the reconstructed MAAT curve shows that the sediments from CO14 are not older than the MIS 7 (191 ka) (Supplementary Fig. S2), deduced by the absence of a contrasting temperature change at the bottom part of the record.

After this validation, we also fitted the MAAT curve to the benthic δ18O LR04 stack[38] using three changes in temperature as tie points (Supplementary Fig. S2, Table S3).The main objective was to refine our chronology for the part beyond the radiocarbon dated section. The selected tie points correspond to the transition from MIS 6 to 5, the onset of the MIS 5e, and the transition from MIS 5/4. Each tie point contains an arbitrary uncertainty of ± 5 ka, acknowledging the differences between both curves.

The age-depth model from the CO3 record was re-calculated based on the five published radiocarbon dates using a Bayesian approach (Supplementary Fig. S4, Table S4) (RBacon[36]). To be able to compare the new proxies and the pollen section from CO14, the CO3 record was aligned to CO14 based on both independent 14C age-depth models (CO3 and CO14). The section with no radiometric control (182 to 752 cm) in the record CO3 was aligned to the CO14 record assuming that both records could have the same sedimentation rate, as they are located close one to another (5 m). Although the two records could hypothetically have similar sedimentation rates due to their proximity, and as shown by the radiocarbon dated section, we acknowledge higher uncertainty at the bottom part of both records, particularly within the CO3 record. For palynological purposes, the two records were merged based on the similarities of the pollen percentages in each data set (CO3 and CO14) (See section on pollen analyses and Supplementary Fig. S6).

XRF analyses

X-Ray Fluorescence (XRF) scans were performed to determine the variability of Al, Si, S, K, Ca, Ti, Fe, Rb and Zr using an ITRAX scanner (Cox Analytical Systems) at CEREGE (France). To optimise the count response of light and heavy elements, scans were performed using a Cr tube set at 30 kV and 40 mA, and a Mo tube at 30 kV and 45 mA, respectively. The spatial resolution of XRF measurements was 1 mm with a counting time of 15 seconds. The ratio of incoherent to coherent responses (inc/coh) obtained with the Mo tube corresponds to the Compton to Rayleigh scattering ratio and represents the relative abundance of light and heavy elements[68,69]. Transmitted X-ray images were obtained using the Mo tube set at 40 kV, 35 mA, the counting time was 500 ms for a resolution of 200 μm.

brGDGTs temperature curve

Mean annual air temperature (MAAT) is based on the analysis of branched glycerol dialkyl glycerol tetraethers (brGDGTs) following the method described by Hopmans et al.[70] and adapted at CEREGE by Davtian et al.[71,72]. A total of 92 sub-samples from CO14, each weighing between 0.1 to 0.9 grams were freeze-dried, ground and extracted using an accelerated solvent extraction system (ASE350 Dionex-Thermo Fisher) at 120 °C, 100 bars and a mixture of dichloromethane and methanol (9 :1, v/v). Extracts were then separated into apolar and polar fractions following the method described by Sanchi et al.[73]. The polar fraction with the brGDGTs was analysed by high-performance liquid chromatography-atmospheric pressure chemical ionisation-mass spectrometry (HPLC-APCI-MS; Agilent 1250 Infinity HPLC coupled with an Agilent 6120 quadrupole mass spectrometer).

We applied the global peat-specific temperature calibration protocol developed by Naafs et al[35]. (MAAT°C = −23.05 + 52.18*MBT′5ME), which led to MAAT values of about 24 °C for samples from the top of the core. This value is higher than the modern MAAT (20 °C), with a seasonal range between 15 and 25 °C at Colônia (Supplementary Fig. S7). However, the brGDGTs-based temperature is still compatible with the modern MAAT, considering the calibration uncertainty (root mean square error (RMSE) ≈ 4.7 °C). The brGDGTs-based temperature values correspond to the summer temperature which is also the wettest period of the year (Fig. 1); both parameters being favourable to the brGDGT producers. Using a recent brGDGTs soil calibration[74] (MAAT°C = −8.57 + 31.45*MBT′5ME; RMSE ≈ 4.8) results in very similar MAAT estimates for the glacial period. By contrast, the application of a calibration for aquatic production of brGDGTs in lakes (e.g. Russell et al.[75]; MAAT = −1.21 + 32.42*MBT′5ME; RMSE ≈ 2.44) led to very high palaeotemperatures values ranging between 27 and 29 °C during MIS 5 and 29 °C for the core top, much higher than the modern MAAT and seasonal range.

In order to estimate the origin of tetraethers at Colônia and selection of the most appropriate calibration curve, we used the ratio of specific GDGTs and brGDGTs as qualitative indicators of terrestrial organic matter versus aquatic production (Supplementary Fig. 7). First, we applied the Branched and Isoprenoid Tetraether (BIT) index, which varies between 0 and 1 (1 being an entirely terrestrial origin)[71,76-78] and the brGDGTs ratio ∑IIIa/∑IIa (∑IIIa = sum of isomers of the ion trace [M + H] + 1050 and ∑IIa = sum of isomers of the ion trace [M + H] + 1036) with values <0.5 considered as terrestrial and >1.2 considered of aquatic origin[79,80]. For both, the lacustrine and peat sections, high BIT index values (0.95 to 1) and very low ∑IIIa/∑IIa ratios (0.04 to 0.22) are the signature of a terrestrial origin for brGDGTs in the Colônia record. In the peat section, brGDGTs production is thus considered as in-situ[35,77], whereas in the lake section accumulated during humid glacial periods, the brGDGTs probably came from the erosion of soils on the borders of the lake or its catchment. This origin is supported by increased PC1 of XRF elemental data indicating siliciclastic sediments deposition during the lacustrine phase (Fig. 2).

Consequently, we applied the peat calibration of Naafs et al.[35] to the entire CO14 core. The calibration leads to cooling of about 8 °C for the MIS 3, MIS 2 and MIS 6. This glacial cooling is somewhat larger than previous estimations based on noble gases in aquifers, pollen in sediments and mountain snow lines in South America for the LGM (ca. 5 °C)[81]. This difference in amplitude may be linked to the uncertainty of the brGDGTs calibration from Naafs et al.[35].(RMSE ≈ 4.7 °C). Using the lake sediment calibration by Russell et al.[75] would lead to lower glacial cooling of about 5–6 °C.

In any case, the overall shape and relative variation of the brGDGTs palaeotemperature curves are similar using different calibrations for peat, soils and lakes[35,74,75]. This implies that the identification of MIS 1, MIS 2, MIS 5 and MIS 6 is reliable, as confirmed by the correlation with changes in other stratigraphic proxies (e.g. XRF elements intensity, charcoal and pollen data). A similar terrestrial origin of the brGDGTs in both, the peat and lacustrine sediments is demonstrated by high values of the BIT index and low values of the ∑IIIa/∑IIa ratio, which justifies the application of the same calibration throughout the whole core, preferably the global peat calibration as it leads to recent temperatures compatible with that of the present day.

Land-sea anomalies

The variations of the land-sea temperature gradient have been calculated by subtracting the temperature anomalies with respect to the core top, for the continental MAAT record (based on brGDGTs), and for the marine SST record (based on Mg/Ca in the core GL-1090)[39]. This index is a semi-quantitative approximation which precision depends on uncertainties of the calibrations used for both proxies (±1.2 and 4.7 °C, respectively) and consequently, only systematic and sustained changes are meaningful and to be considered in our study (e.g. contrasts between MIS 5, MIS 2, 3, 4 and MIS 6) (Supplementary Fig. S8).

Pollen analyses

Pollen analyses were based on two parallel cores separated by a distance of 5 m, CO3[17] and CO14. Based on the previously described core alignment, and acknowledging the uncertainties from the bottom part of the record CO3, we considered the matching point to be 105 ka (CO3:726 cm and CO14: 690 cm) based on arboreal pollen, Poaceae and Araucaria percentage matching curves. We have excluded only 6 samples in the CO3 data set for graphical purposes (Supplementary Fig. S6). Consequently, pollen samples from the interval 180 to 105 ka belong to the CO14 record and pollen samples from 105 ka towards the present belong to the CO3 record[17].

A total of 57 subsamples were analysed for the CO14 record. We applied an uneven subsampling approach based on knowledge of eight test samples (1263, 1277, 1293, 1299, 1301, 1313, 1327 and 1337 cm depth) of 1 cm3 in volume. The samples were collected from the clayey section of the record with high pollen content; the remaining unanalysed section was subsampled using a smaller volume, 0.5 cm3. Each sample was treated using standard methods for pollen analysis[82]. Pollen was counted to a minimum of 300 terrestrial pollen types.

A percentage peak of Podocarpus from the record CO3 centred at 83 ka (lasted ca.1 ka) is associated to specific climate features with a sharp decrease in MAAT, land-sea gradient and an increase in pCO2. Pollen assemblages are composed of 9–16% of Ilex, 2–8% of Asteraceae and 2–9% of Myrtaceae, and do not relate to any modern analogue of the Atlantic forest today[34]. This suggests that Podocarpus expanded locally on the lake shore and therefore its peak is only highlighted with an asterisk in Figs. 2 and 3.

All illustrations and calculations were done using the Psimpoll program[83] and R[56].

Charcoal analyses

A total of 147 subsamples were sampled for microcharcoal analyses according to the methodology of Daniau et al.[84] using a series of chemical laboratory treatments combining hydrochloric acid (HCl), nitric acid (HNO3), hydrogen peroxide (H2O2) and hydrofluoric acid (HF) to remove calcium carbonates, organic matter and siliceous material. The peat section (850 to 5 cm) was treated differently due to the very high concentration of organic remains. A total of 0.5 g of dry peat was sieved with water. Next, 0.2 g of dry sediment (<150 µm) per sample was chemically treated, excluding the step using HF. Microcharcoal particles were counted manually at x 500 magnification.

Supplementary Material.