Volcanism and climate change as drivers in Holocene depositional dynamic of Laguna del Maule (Andes of central Chile – 36° S)

Late Quaternary volcanic basins are active landscapes from which detailed archives of past climate and seismic and volcanic activity can be obtained. A multidisciplinary study performed on a transect of sediment cores was used to reconstruct the depositional evolution of the high-elevation Laguna del Maule (LdM) (36^∘ S, 2180 m a.s.l., Chilean Andes). The recovered 5 m composite sediment sequence includes two thick turbidite units (LT1 and LT2) and numerous tephra layers (23 ash and 6 lapilli). We produced an age model based on nine new ¹⁴C AMS dates, existing ²¹⁰Pb and ¹³⁷Cs data, and the Quizapú ash horizon (1932 CE). According to this age model, the relatively drier Early Holocene was followed by a phase of increased productivity during the mid-Holocene and higher lake levels after 4.0 ka cal BP. Major hydroclimate transitions occurred at ca. 11, 8.0, 4.0 and 0.5 ka cal BP. Decreased summer insolation and winter precipitation due to a southward shift in the southern westerly winds and a strengthened Pacific Subtropical High could explain Early Holocene lower lake levels. Increased biological productivity during the mid-Holocene (∼8.0 to 6.0 ka cal BP) is coeval with a warm–dry phase described for much of southern South America. Periods of higher lake productivity are synchronous to a higher frequency of volcanic events. During the Late Holocene, the tephra layers show compositional changes suggesting a transition from silica-rich to silica-poor magmas at around 4.0 ka cal BP. This transition was synchronous with increased variability of sedimentary facies and geochemical proxies, indicating higher lake levels and increased moisture at LdM after 4.0 ka cal BP, most likely caused by the inception of current El Niño–Southern Oscillation and Pacific Decadal Oscillation (ENSO–PDO) dynamics in central Chile.

1 Introduction

The Andes of central Chile (or southern Central Andes; SCA hereafter) are prone to drought, and future scenarios of global warming show them becoming drier as winter droughts are more recurrent (Falvey and Garreaud, 2009; Stocker et al., 2013; Boisier et al., 2016). Many studies have documented major changes in the terrestrial ecosystems and atmospheric and oceanic circulation associated with dry/wet periods during the Holocene on the western slope of the Andes in Chile (e.g., Jenny et al., 2002; Lamy et al., 2001; Martel-Cea et al., 2016; Latorre et al., 2007, 2006; Frugone-Álvarez et al., 2017; Fletcher and Moreno, 2011; Kim et al., 2002; Kaiser et al., 2008), but the nature, regional distribution and timing of these variations in the Mediterranean-type regions of central Chile are not yet fully understood. Several coupled atmosphere–ocean mechanisms, including variations in tropical Pacific sea surface temperature (SST) gradients, poleward displacement of the southern westerly wind (SWW) and changes in the intensity and position of the South Pacific Subtropical High (SPSH), have been hypothesized to explain this climate variability (Jenny et al., 2002; Lamy et al., 2010; Valero-Garcés et al., 2005; Rein et al., 2005; Liu, 2014; Haug et al., 2001). Furthermore, little is known about the long-term evolution of large-scale climate change in these regions of the Southern Hemisphere, owing to the lack of high-resolution paleoclimatic datasets (Villalba et al., 2009).

The SCA are also one of the most active volcanic and seismic zones on Earth due to the subduction of the Nazca plate underneath the South American plate (Dewey and Lamb, 1992). Recent volcanic activity, with one major eruptive event in the last 200 a from the Quizapú Volcano (Plinian eruption with a volcanic explosivity index (VEI) =5; Fontijn et al., 2014), underlines that this is a regional hazard for the society and the economy of central Chile (Hildreth and Drake, 1992). With ongoing uplift >60 m of permanent deformation since the last deglaciation, the magmatism in Laguna del Maule (LdM) has been reconstructed based on geochemistry and the dating of surface volcanic formations (Hildreth et al., 2010; Andersen et al., 2017; Singer et al., 2018). In recent decades, the LdM volcanic field has shown extremely high rates of deformation (Andersen et al., 2012; Feigl et al., 2014; Singer et al., 2018; Wespestad et al., 2019) that can be interpreted as a response to shallow magma intrusion, thus indicating intensified volcanic activity and hazard (Fig. 1).

Figure 1 Regional context of LdM in the southern Central Andes of central Chile. (a) Schematic map of the principal low-level atmospheric flows over South America. SWW – southern westerly wind; SPA – southeast Pacific anticyclone; TW – trade winds; ITCZ – Intertropical Convergence Zone; LLJ – subtropical low-level jet stream; SACZ – South Atlantic Convergence Zone; and ACC – Antarctic Circumpolar Current. Light green (tropical rainfall), orange (stratocumulus and cold SST) and light blue (midlatitude precipitation) areas represent the influence of atmospheric precipitation associated with each forcing, modified from Garreaud et al. (2009). The maps (b, c) show the main locations of the paleoclimate records discussed in the text with a digital elevation model (STRM30; Becker et al., 2009; Smith and Sandwell, 1997). White circles are lacustrine records: 2 – El Junco lake (Zhang et al., 2014); 3 – L. Pallcacocha (Moy et al., 2002); 6 – L. Queshquecocha (Stansell et al., 2013); 7 – L. Pumacocha (Bird et al., 2011); 9 – L. Titicaca (Fornace et al., 2014); 10 – L. Chungará (Moreno et al., 2007); 11 – L. Miscanti (Valero-Garcés et al., 1996, 1999); 12 – L. Negro Francisco (Grosjean et al., 1997; Valero-Garcés et al., 1999); 17 – L. Chepical (de Jong et al., 2013); 20 – L. Matanzas (Villa-Martínez, 2002); 21 – L. Aculeo (Jenny et al., 2002, 2003); 22 – L. Vichuquén (Frugone-Álvarez et al., 2017); 23 – L. Tagua Tagua (Valero-Garcés et al., 2005); 25 – L. Laja (Urrutia et al., 2010); 26 – Lanalhue and Lleu Lleu lakes (Stefer et al., 2010); and 27 – Pichilafquén lake (Jara and Moreno, 2014). Green pentagons represent coastal peat swamp records: 14 – Ñague and Quereo (Maldonado and Villagrán, 2002; Villagrán and Varela, 1990); 15 – Palo Colorado (Maldonado and Villagrán, 2006); and 19 – Quintero ((Villa-Martínez and Villagrán, 1997). Blue diamonds represent marine records: 1 – Cariaco Basin ODP1002 core (Haug et al., 2001); 8 – SO147-106KL (Rein et al., 2005); 13 – GEOB7139-2 (De Pol‐Holz et al., 2007); 16 – ODP 1233 and GeoB3313-1 cores (Muratli et al., 2010a, b; Lamy et al., 2002); 18 – GEOB3302-1 and GIK17748-2 cores (Kim et al., 2002); and 24 – ODP 1234 and ODP 1235 cores (Muratli et al., 2010a, b). Yellow squares represent speleothem records: 4 – Santiago cave (Mosblech et al., 2012); and 5 – Tigre Perdido cave (van Breukelen et al., 2008). The software used to build the map was QGIS 2.8.

Here, we present a seismic survey and a multiproxy sediment core transect study from LdM to investigate depositional variability in a lake located in one of the most active volcanic fields in the world. Although intense volcanic and seismic activity poses a challenge to paleoenvironmental reconstructions, they also provide an opportunity to investigate the complex interplay of these factors in the lake evolution. Geochemical proxies are interpreted in terms of endogenic productivity and detrital input variation associated with climate variability and the influence of volcanic activity in the lake basin. Pollen proxies provide information about regional and local vegetation dynamics. The LdM lacustrine sedimentary sequence contains a detailed archive of volcanic history and exemplifies how volcanic activity can impact lake processes, not only in terms of tephra input but also on local geochemical and biological cycles. As the sequence spans the entire Holocene, the LdM record also provides one of the few opportunities to test previous hypotheses regarding the coupled atmosphere–ocean mechanisms that controlled the climate evolution of central Chile during the Holocene (Jenny et al., 2002; Lamy et al., 2001; Villagrán and Varela, 1990; Valero-Garcés et al., 2005; Rutllant and Fuenzalida, 1991; von Gunten et al., 2009), especially with respect to the relative contributions of the SWW, SPSH, volcanic activity and tropical Pacific SST gradients during known rapid climate change events (Mayewski et al., 2004).

2 Regional setting

The LdM volcanic field is part of the Southern Volcanic Zone (SVZ, 39–34^∘ S) in the SCA (36^∘ S–70^∘30^′ W, 2180 m a.s.l.; Fig. 1c), which is one of the most active volcanic and seismic zones on Earth (Feigl et al., 2014). The most voluminous Quaternary eruption of the volcanic field occurred ca. 950 ka and formed part of the LdM basin as a consequence of the collapse of an 80 km² elliptical caldera (Bobadilla caldera). Between 336 and 38 ka, small rhyolitic eruptions took place, culminating in a ring of 36 postglacial rhyodacite and rhyolite coulees and domes over the last 25 kyr (Andersen et al., 2012; Feigl et al., 2014; Singer et al., 2018, Fig. 2).

Figure 2 Bathymetric and geological map of LdM volcanic complex (based on Hildreth et al., 2010; Andersen et al., 2017; Singer et al., 2018). White, green, gray and red circles indicate coring sites, volcanic vents and sample locations for ⁴⁰Ar,³⁹Ar and ³⁶Cl dating (Andersen et al., 2017) and El Piojo lake, respectively. The cores used to determine the reservoir effect correction are derived from site 3 (LEM13-3D-G: ¹⁴C age for the wood) and site 12 (LEM11-3A-G: ¹⁴C age for macrophytes). Gray circles showing eruptions (>25 kyr) composed principally of the ignimbrite of Cajones de Bobadilla (igcb; ∼950 kyr), rhyolite of Cerro Negro (rcn; ∼447 kyr), rhyodacite of Domo del Maule (rddm; ∼114 kyr), basalt of El Candado (bec; ∼63 kyr) and the andesite of Arroyo Los Mellicos (aam; ∼26 kyr). Eruptions (<25 kyr) include the rhyolite east of Presa Laguna del Maule (rep; 25.7 kyr), rhyolite of Loma de Los Espejos (rle; ∼19 kyr), rhyolite of Cari Launa (rcl; <3.3 kyr), rhyolite south of Laguna Cari Launa (rsl; 3.3±1.2 kyr), rhyolite of Arroyo de Sepúlveda (ras; 20–19 kyr), rhyolite of Cerro Barrancas (rcb; multiple flows 11.4–1.9 kyr), rhyolite of Colada Divisoria (rcd; 2.1±1.3 kyr), rhyolite of Colada Las Nieblas (rln; Late Holocene), rhyodacite of Arroyo de la Calle (rdac; $\sim \mathrm{20}\pm \mathrm{1.2}$ kyr), rhyodacite of Colada Dendriforme (rdcd; 8±0.8 kyr), rhyodacite of the northwest coulee (rdcn; 3.5±2.3 kyr), rhyodacite of Laguna Sin Puerto (rdsp; <3.5 kyr), rhyodacite west of Presa Laguna del Maule (rdop), andesite of Laguna Sin Puerto (asp; <3.5 kyr) and the younger andesite of the west peninsula (apj; 21±3.4 kyr). The age uncertainties are 2σ.

The basin (∼300 km²) is irregularly shaped (Fig. 2) due to both tectonic/volcanic processes and glacial erosion (Singer et al., 2000). The maximum elevation of the watershed reaches 3940 m a.s.l. at Cerro Campanario and the mean elevation is 2450 m a.s.l. (Fig. 2). The lake has ∼50 m maximum water depth, and the Maule River is its only outlet (Fig. 2). The lake level is mainly controlled by the influx of snowmelt through small inlets and runoff. A dam was completed at the lake outlet in 1957, and, as a result, water storage increased from 600 to ca. 2000 hm³ with a 30 m lake level increase (Carrevedo et al., 2015). Despite limited available limnological data, the local climate and basin morphology suggest that LdM is a dimictic lake (Frugone-Álvarez, 2016). The temperature profile in March 2013 (austral autumn) revealed that LdM was slightly stratified, with a gradual temperature decrease from 12 ^∘C at the surface to 10.5 ^∘C near the bottom (Fig. S3b in the Supplement). The waters were oligotrophic with relatively high alkalinity ([CO₃Ca]=424 mg kg⁻¹), pH between 7.0 and 8.4, and a low concentration of salts with nitrate and sulfate values reaching 50 and 7000 µg L⁻¹, respectively (Frugone-Álvarez, 2016). Relatively high alkalinity could be related to the dissolution/hydrolysis of volcanic rocks (Singer et al., 2014). Total phosphate was very low, 10 µg L⁻¹, which may potentially limit biological activity. The δD_[VSMOW] and δ¹⁸O_[VSMOW] values from small inlets and water samples at different depths in the lake suggest that lake waters are derived mainly from precipitation and snowmelt with little deviation from the local meteoric water line (Fig. S3b).

Regional climate is characterized by ∼1700 mm annual rainfall and a large seasonal temperature range (Fig. S4). Winter precipitation, mostly as snow, is associated with the incursion of cold fronts that migrate with the SWW (Garreaud, 1992; Falvey and Garreaud, 2007). Summers are dry due to the strong influence of the SPSH, and easterly storms are infrequent (Viale and Garreaud, 2013). On interannual to interdecadal timescales, rainfall in subtropical central Chile is sensitive to variations in the El Niño–Southern Oscillation (ENSO), the Southern Annular Mode (SAM) and the Pacific Decadal Oscillation (PDO) (Jacques-Coper and Garreaud, 2014; Garreaud et al., 2009; Montecinos and Aceituno, 2003; Pittock, 1980; Quintana and Aceituno, 2012). These climate modes have a strong influence on the patterns of variability of the snowpack and streamflow in central Chile along with the occurrence of extreme events, such as protracted drought or torrential rains (Masiokas et al., 2010, 2012).

Modern vegetation surrounding LdM is dominated by sub-shrubs and cushion species corresponding to the high Andean shrubland belt – between 2000 and 2500 m a.s.l. (Luebert and Pliscoff, 2006). At higher elevations (>2500 m a.s.l.), the high Andean steppe is characterized by scattered and/or scarce vegetation dominated by herbs and cushion species like Oxalis adenophylla and Pozoa coriacea. Low Andean shrublands occur at elevations between 1500 and 2000 m a.s.l. and include Chuquiraga oppositifolia, Discaria articulata, Laretia acaulis, Berberis empetrifolia and Discaria chacaye. Deciduous forest is found below 1500 m a.s.l. (Fig. S5).

3 Methods

3.1 Seismic surveys and coring

The LdM sedimentary basin was surveyed with an EdgeTech SB-424 sub-bottom multifrequency profiler using a frequency range of 2–10 kHz. Approximately 20 km of seismic reflection data were acquired in the northern areas of the basin with a denser grid network close to the coring sites. The data were processed and interpreted using the EdgeTech Discover SB 3200-XS software. The cores were linked to the seismic data with a depth–time conversion, assuming an acoustic velocity of 1500 m s⁻¹ and using the density values and the P-wave velocity measured in the cores with a Geotek Multi-Sensor Core Logger (MSCL). Coring campaigns were organized during the summers of 2011, 2012 and 2013. We recovered 17 short cores and over 40 m of long sediment cores at several sites using a hammer-modified Uwitec gravity corer and a Uwitec platform with a percussion piston corer (Fig. 2; Table S1 in the Supplement).

3.2 Core analyses

Sediment cores were transported to the IACT-CSIC laboratory in Granada, Spain, to measure physical properties (P-wave velocity, magnetic susceptibility, electric resistivity and gamma density) with a 1 cm resolution by a Geotek MSCL. The cores were then imaged with a Geotek Single Track Core Imaging System (MSCL-CIS) to 50 µ pixel size at the Pyrenean Institute of Ecology (IPE-CSIC) labs, Spain.

Sedimentary facies were defined and characterized based on macroscopic (color, bedding features, grain size, lithology and sedimentary textures) and microscopic (smear slides) criteria formulated by Schnurrenberger et al. (2003). A composite sedimentary sequence (4.6 m long) was constructed with the short and long cores at site 3, LEM13-3A, LEM13-3E and LEM13-3B, and also LEM13-2C at site 2, which was used for geochemical and pollen analyses (Table S1).

Elemental geochemistry included total carbon (TC), total inorganic carbon (TIC), total organic carbon (TOC = TC − TIC) and total sulfur (TS) at 1 cm resolution performed with a LECO elemental analyzer and total nitrogen (TN) at 5 cm resolution in a vario MAX CN elemental analyzer at the IPE-CSIC. Biogenic opal (BioSi) was measured using the wet alkaline leaching procedure, and dissolved silica was photometrically determined using the continuous flow analyzer AutoAnalyzer Technicon II at the Instituto de Investigaciones Marinas (IIM-CSIC), Spain (Bernárdez et al., 2005). Stable carbon and nitrogen isotope analyses were performed on bulk organic matter (OM) from sediment samples of 20–30 mg (dry weight) at 5 cm intervals and measured in the Laboratory of Biogeochemistry and Applied Stable Isotopes (LABASI) of the Pontificia Universidad Católica de Chile using a Thermo Delta V Advantage IRMS coupled with a Flash2000 elemental analyzer.

An AVAATECH X-Ray Fluorescence II core scanner with an Rh X-ray tube from the University of Barcelona was used to obtain X-ray fluorescence (XRF) logs from the LEM11-3A, LEM11-3E and LEM11-3B cores at 4 mm resolution. Only elements with mean values higher than more than 1000 counts per second (cps) were used in the statistical analyses. We use a robust principal component analysis for dimensionality reduction (Fig. S8). Three element ratios have been calculated to describe changes in the redox conditions (Fe∕Mn), carbonate (Ca∕Ti) and organic productivity (Br∕Ti) (Carrevedo et al., 2015; Frugone-Álvarez et al., 2017; Naeher et al., 2013; Moreno et al., 2007; Fuentealba et al., 2020). Thirty-one discrete samples were measured using inductively coupled plasma–optical emission spectrometry (ICP-OES) at the CEBAS-CSIC laboratory in order to analyze the relationship between quantitative elements and XRF semi-qualitative measurements (Fig. S10). Sixty-four pollen samples were prepared according to standard techniques (Faegri et al., 1989), identified and counted using reference collections available at the CEAZA, Chile (Carrevedo et al., 2015). The pollen ratio was calculated using the formula $(a-b)/(a+b)$, where a is Poaceae abundances (%) and b is Ephedra abundances (%) (Maher, 1963, 1972). In this way, values near 1 in the pollen ratio suggest an increase in humidity, while values near −1 suggest drier conditions. The data were analyzed and plotted using QGIS software 2.8 (QGIS Development Team, 2009) and R language environment 3.5.1 (R Core Team, 2018), using the packages analogue (Simpson, 2007), dplyr (Wickham et al., 2018), vegan (Oksanen et al., 2018), ggplot2 (Wickham, 2016) and pcaPP (Filzmoser et al., 2018).

3.3 Chronology

The chronology for the LdM sequence was obtained by ²¹⁰Pb and ¹³⁷Cs dating (Carrevedo et al., 2015) and nine new accelerator mass spectrometry (AMS) radiocarbon dates on wood, terrestrial macrofossils, algal macrofossils and bulk sediment (Table 1). Radiocarbon ages were determined by AMS ¹⁴C at the Poznan Radiocarbon Laboratory, Poland, the DirectAMS laboratories and UC Irvine (Keck radiocarbon facility), USA. Waters from the hypolimnion and modern living littoral macrophytes were also sampled to assess modern ¹⁴C reservoir effects (Table 1). We used the R package Bacon v2.2 (Blaauw and Christen, 2011) for age–depth modeling and to establish the deposition rates along the core after removal of instantaneous depositional events such as tephra layers and turbidites (∼100 cm from the original sequences was removed to calculate the age model). Radiocarbon data are reported as radiocarbon age in years before present (relative to 1950 CE). Radiocarbon ages were calibrated using the Southern Hemisphere calibration curve (SHCal13) applying the reservoir effect in the ages of aquatic organic matter inside the R package Bacon (Hogg et al., 2013). The short core LEM12-3B was sampled in the field for ²¹⁰Pb∕¹³⁷Cs dating every 0.5 cm for the uppermost 20 cm and at 1 cm for the lower section (Table S3). Analyses were performed at the St. Croix Watershed Research Station Laboratory, Science Museum of Minnesota, USA. Dates and sediment accumulation rates were established using the constant rate of supply (CRS) model (Appleby and Oldfield, 1978). According to the stratigraphic correlations between short cores and the ²¹⁰Pb∕¹³⁷Cs age model, we assumed that T1 is the 1932 Plinian eruption (Fontijn et al., 2014). The assumed Quizapú ash layer (1932 CE) and the appearance of Pinus pollen served as further chronological markers, as Pinus radiata was not forested in large-scale plantations in Chile until 1931 CE (Ley 4.363, 30 June 1931, Santiago, Chile; Table S3). Dates for the LdM eruptive history based on ¹⁴C dating of paleosols, ⁴⁰Ar∕³⁹Ar dating of lava flows and ³⁶Cl dates or paleoshoreline exposures were compiled from the literature (Andersen et al., 2017) to better constrain the timing of the depositional evolution of LdM.

Table 1 Radiocarbon dating of LdM sequence. We used only the macrophyte samples for the construction of age–depth model with a reservoir correction of 4700 a. The table shows the median age with reservoir effect corrected from the Bayesian age–depth model and uncertainties (lower and upper error).

^* Indicates the AMS radiocarbon dates without ¹⁴C reservoir effect. ${}^{**}$ Dating excluded from the age model. ${}^{***}$ Radiocarbon age of dissolved inorganic carbon (DIC) from hypolimnion to 20 m depth and modern living littoral macrophytes used to assess modern ¹⁴C reservoir effect.

4 Results and interpretation

4.1 Sedimentary facies

Sediments in LdM are composed of diatom-rich muds, silts and oozes with abundant interspersed volcanic facies (cryptotephra and tephra like ash and lapilli). The presence of centimeter-thick macrophyte and diatom-rich silt facies characterizes deposition in littoral, shallower (<25 m water depth) cores (e.g., LEM11-3A, ∼24 m water depth) (Figs. 2 and 3). In more distal, deeper areas (>25 m water depth; sites 1, 2 and 3), sediments are made up of banded to laminated organic and diatomaceous silts and include several thick massive layers interpreted as lacustrine turbidites. LdM cores from sites 1, 2 and 3 have been correlated using sedimentary facies, volcanic layers and physical properties (magnetic susceptibility (MS) and density) (Fig. 3). Six distal hemipelagic facies, 2 mass wasting depositional events and 29 volcanic layers (6 lapilli and 23 ash) have been identified in the sedimentary sequences of LdM (Figs. 3 and 4).

Figure 3 West–east LdM core correlation in the northern areas based on lithostratigraphic and sedimentological criteria. White circles indicate coring site. Magnetic susceptibility (MS; $\times {\mathrm{10}}^{-}\mathrm{8}\left({\mathrm{m}}^{\mathrm{3}}{\mathrm{kg}}^{-\mathrm{1}}\right)$), total organic carbon (TOC, %) and the occurrence of tephra layers. The stratigraphic units are marked in different colors: unit 1 – yellow; unit 2 – light brown; unit 3 – brown; unit 4 (LT1) – light red; unit 5 – light yellow; and unit 6 (LT2) – red.

Figure 4 Sedimentary facies and sedimentological units in the LdM sequence for site 3. Six lithostratigraphic units (U) and three main facies groups: lacustrine facies (D1 to D6); lacustrine turbidites (LT1 and LT2); and volcanic facies (lapilli L1 to L5 and tephras T1 to T18). The lacustrine facies have been classified according to elemental composition: magnetic susceptibility (MS; $\times {\mathrm{10}}^{-}\mathrm{8}\left({\mathrm{m}}^{\mathrm{3}}{\mathrm{kg}}^{-\mathrm{1}}\right)$), percentage of the total organic carbon (TOC), total sulfur (TS), total inorganic carbon (TIC), biogenic silica (BioSi), atomic TOC∕TN values, δ¹⁵N and δ¹³C values in per mill (standardized with N2-Air and Vienna Pee Dee Belemnite (VPDB), respectively), and XRF ratios as proxies for redox conditions (Fe∕Mn) and organic (Br∕Ti) productivity.

The upper five units in site 1 are similar to site 3 (Fig. 3), suggesting that the composite stratigraphy of site 3 is representative of the sedimentary infill in the northern areas of the basin. Core LEM13-1B (31 m water depth) reached 5.8 m beneath the lake floor, recovering 2.5 m of lacustrine sediments and volcanic facies below laminated unit 5. Unit 6 (LT2) in the shallow, western site 1 is composed of gravels and sands with an abundant silty matrix, topped by a thin, homogeneous layer. Sediments below unit 6 in site 1 have been grouped in stratigraphic unit 7 (Fig. 3), including coarse clastic facies, fine lacustrine laminated sediments, thin lapilli (L6, 25 cm) and five ash layers (T19 to T23). The base of the core in site 1 is composed of coarse breccia with angular centimeter-long volcanic clasts in a volcanic matrix that likely represents brecciated volcanic facies emplaced in the distal areas of the lake. Laminated diatomaceous facies cover this basal breccia with an intercalated thick, black, mafic-rich ash (T23). The presence of homogeneous coarse sands with a clay–silt matrix and coarse gravels with an unsorted sandy matrix and rounded clasts suggests fluvial/alluvial depositional processes in the lake.

4.1.1 Hemipelagic diatomaceous facies

Lacustrine diatomaceous facies have been classified according to elemental composition (TOC, TS, BioSi, TIC), grain size (textures), sedimentary structure (lamination) and biological composition (diatoms and macrophytes) in three groups: (i) banded silts with low TOC (<2 %) (D1 and D2), (ii) banded organic-rich silts (>2 % TOC) (D3, D4 and D5) and (iii) laminated organic-rich silts (D6) (Figs. 4 and S6; Table 2). Banded, organic-rich silts occur in the upper half of the sedimentary sequence at sites 1, 2 and 3 (facies D1 to D5), whereas laminated facies D6 occur in the lower part.

Table 2 Sedimentological and compositional characteristics and depositional environment of LdM lacustrine facies.

Table 2 summarizes the main characteristics of the LdM lacustrine facies. The finer grain size of facies D1 and D2 and the absence of littoral components (i.e macrophyte remains) indicate deposition in relatively deep water. Coarser grain size and the abundance of macrophyte remains suggest a more littoral depositional area for facies D3 compared to D1 and D2. Facies D3, D4 and D5 are organized in decimeter-thick sequences and are macrophyte-dominated (D5, D3) to diatom-dominated environments (D4) (Fig. S6). Laminated facies D6 have the highest TOC and TS values up to 5.5 % and 7.0 %, respectively, and contain euhedral crystals of endogenic calcite about 10 µm in size.

Biogenic silica concentrations are between 5 % and 26 % with a mean of 20.5 %. A boxplot of BioSi distribution in the sedimentary facies (Figs. 4 and S6h) shows the highest values in banded facies D3, D4 and D5 and relatively lower values in laminated (D6) and core top sediments (D1 and D2). Well-defined low BioSi values occur at the base (432–433 cm), middle (364–366 cm) and top (313–316 cm) of unit 5, some of them associated with volcanic facies (Fig. 4).

4.1.2 Mass wasting deposits and lacustrine turbidites

Two thick intervals of homogeneous sediments occur in the middle (Facies LT1) and at the base (Facies LT2) of the sequence (Fig. 3 and Table 3). They show relatively constant XRF intensities and δ¹⁵N values but with different δ¹³C and TOC∕TN value distributions between LT1 and LT2 (Figs. 4, S6, S9, S10 and S12a, b), pointing to a mixture of different terrestrial, macrophyte and phytoplankton organic matter (Meyers, 2003, 1994) or metallic CO₂ input that is less negative than atmospheric CO₂. LT1 presents a banded, coarse-grained, 30 cm thick basal interval and a thicker (up to 1 m) homogeneous upper interval (Figs. 4, S6 and S9). LT2 is composed of finer black silts than LT1, while the latter is browner and coarser and has more abundant macrophyte remains and millimeter-size pumice clasts than the former. Although the core at site 3 did not reach the base of this unit (Fig. 3), the lower part of LT2 shows a small increase in MS (Figs. 4 and S6), which is suggestive of graded texture. The upper limits for both LT1 and LT2 are sharp. LT2 is overlain by finely laminated, diatom-rich facies and LT1 by a tephra layer (T7).

Table 3 Depositional and compositional characteristics of LdM mass wasting facies.

The sedimentological features (gradational textures with basal layering and homogeneous at the top) and the homogeneous chemical composition are key criteria for establishing these layers as mass wasting deposits emplaced in the deeper parts of LdM. Seismic profiles (Fig. 6) show some massive deposits – sublacustrine mass-transport deposits (MTDs), landslides – related to LT1 and focused on the deepest parts of the basin, while the lacustrine facies have a much more draped appearance (Figs. 6b, S1 and S2), supporting this interpretation.

This double structure within the turbidites (banded and massive in LT1 and high basal MS in LT2) has been described in lacustrine turbidites in southern central Chile (Van Daele et al., 2015; Moernaut et al., 2014, 2009), the Alps (Lauterbach et al., 2012) and the Pyrenees (Corella et al., 2014). Additionally, the higher values of TIC, BioSi, δ¹³C and the abundance of macrophyte remains in LT1 point to a higher contribution of littoral sediments (Figs. 4, S6 and S12). Deposition of LT1 occurred after the complex volcanic event centered around Lapilli 3 (T9–L3–T8; Figs. 3 and 4; see below), and that could explain the higher content in pumice clasts. The association with volcanic layers at the base and top of LT1 suggests that volcanic activity and related local seismicity could have been the main trigger for sediment destabilization on the lake margins. In this context of increased seismicity and volcanism, the alteration of pyroclasts (mainly glass/pumice) favors a higher alkalinity that increases the higher carbonate precipitation. Possible deep CO₂ contributions would translate into less negative δ¹³C values in primary production (Fig. 4).

4.1.3 Volcanic facies

At site 3, distinctive volcanic facies occur as lapilli (5L layers) and ash (18T layers) (Table 4). Ash layers have been classified into two main groups (Fig. S6) according to macroscopic features (grain size, color, depositional structures), microscopic observation of smear slides and compositional data from XRF analyses:

• 1.

  Dark colors (T1, T2, T4, T5, T11 and T14) represent those with relatively high content of mafic minerals and low quartz and plagioclase (Frugone-Álvarez, 2016). Semi-quantitative XRF data show higher Sr, Ti and Fe and lower Si (and thus lower Si∕Ti ratio), K and Rb (Figs. S6a, b, S8a). T4 has the highest Ca values. This group includes most layers in the upper part of the sequence after Lapilli 1 (T1 to T5).

• 2.

  Light colors (gray) (T3, T6, T7, T8, T9, T10, T12, T13, T15, T16, T17 and T18) represent those with relatively high glass content and Si (high Si∕Ti ratio), low Fe and Sr, and variable Rb (Figs. S6a, b, S8a). This group mainly includes layers in the lower part of the LdM sequence.

These ashes are laterally continuous throughout the northern and central areas of the basin, where we have a transect of cores. The presence of some depositional textures as normal grading and lamination and the absence of current-derived features are all indicative of fallout deposits from regional and local volcanic activity (Fisher and Schmincke, 1984; Sáez et al., 2007).

Table 4 Depositional and compositional characteristics of LdM volcanic facies.

In the upper part of the sequence, centimeter-thick, volcanic-rich, dark brown to gray homogeneous silts with gradational upper and lower limits occur. They are composed of a mixture of volcanic (glass, silicate crystals, mafic rocks) and organic (diatoms, macrophytes, amorphous) components. They are interpreted as reworked tephra material, deposited in the distal areas of the lake by fluvial and runoff processes, affecting the watershed and washing out the volcanic material after the eruptions.

Lapilli occur as thick (from 6 to 20 cm; see Fig. 3) layers with up to 1 cm long pumice clasts; low matrix content; sharp upper and lower boundaries; and massive, banded or graded (both fining and coarsening upward) textures. Texture, composition and magnetic susceptibility identify two groups: the younger deposits L1 and L2 have more heterogeneous clasts in composition and size, and they have higher Sr, Ca and Fe values (MS more than $\mathrm{200}\times {\mathrm{10}}^{-\mathrm{8}}$ (m³ kg⁻¹)) than the older deposits L3, L4 and L5 (MS less than 100 SI) (Figs. 3, 4, S7c, d and S9). The lapilli layer L3 shows the largest variability in thickness from 2 cm at the eastern sites 3E and 3B to 15 cm at the westernmost site 3A (Fig. 3), suggesting a western source for this volcanic material. It is also the only lapilli layer with a basal ash layer (T9).

4.2 Lithological and seismic stratigraphy

Five seismic units have been defined based on the characteristics of seismic reflectors and correlated with lithostratigraphic units at sites 2 and 3 (Sect. S2.2). Seismic penetration reached up to 10 m in most areas, although only in the upper part can well-defined reflectors be traced all over the northern part of the basin (Fig. 6). The total thickness of the sedimentary infill and morphology of the substratum of LdM could not be established since the substratum–infill boundary could not be seismically imaged. The available seismic profiles show the complex basin structure of LdM, with several sub-basins with variable accumulation rates, a number of faults affecting the lower part of the sediment sequence and lateral changes in thickness and physical properties (Figs. S1 and S2). Using the magnetic susceptibility and density measured in the sediment cores and considering a constant acoustic velocity of 1500 m s⁻¹, the main reflectors in the seismic profiles were correlated with the lithostratigraphic sequence (Fig. 6).

Unit S1 uniformly drapes across the entire northern basin, without significant modification in thickness or physical properties even in the northeastern bay (LEMB10). According to the seismic-to-core correlation of the seismic profile LEMA3 (Fig. 6), unit S1 encompasses U1 to U3 lithostratigraphic units (Fig. 3). Seismic units S2 and S4 have transparent to low-amplitude seismic facies, and these correlate with stratigraphic unit 4 (lacustrine turbidite; LT1) and unit 6 (LT2), respectively. Deposits without a consistent internal structure and with a positive topography occur close to site 3, and they are related to units 4 and 6 (Fig. 6), which can probably be interpreted as mass-transport deposits (MTDs)/landslides. Some profiles show a complex internal structure in marginal areas in S2 (Sect. S2.2). Seismic unit S3 is characterized by closely spaced, well-defined although not very high-amplitude reflections, wavier than in S1. Unit S3 correlates with stratigraphic unit U5.

4.3 X-ray fluorescence (XRF) core scanner geochemistry

The first two principal components (PCs) of the XRF core scanner dataset for lacustrine facies (volcanic facies excluded) explain more than 80 % of the variance (Table S2 and Fig. S8). The eigenvector associated with the highest eigenvalue (5.1) defines two main groups of loadings: (i) Ti (0.33), Rb (0.39), Zr (0.39), K (0.31), Ca (0.35) and Sr (0.32), and (ii) As (−0.30) and Fe (−0.23). The second largest eigenvalue (1.99) is defined mainly by Br (0.50), Mn (0.62) and Ti (−0.40). The first eigenvector shows positive coefficients with higher contents of clastic-related elements, and, consequently, they are interpreted as increasing clastic delivery to the lake (mostly silicates from the volcanic watershed). Negative Z₁-score values (standardized scoring coefficients) are associated with higher values of As and Fe and interpreted as more dominant oxic bottom conditions during periods of lower sediment delivery to the lake (Table S2). This component could also reflect water/sediment input into the system, with positive values of Z₁ score during more humid periods (higher sediment delivery) and negative values during less-humid periods (lower sediment delivery). The Z₂ score is interpreted as an index of bioproductivity with positive coefficients indicative of increased total organic productivity (and/or preservation) and negative values during periods with lower bioproductivity/preservation (Table S2 and Fig. S8).

The Fe∕Mn ratio has the lowest values in laminated facies D6 (lowest Fe, highest Mn; Figs. S9 and S10) and the highest values in banded facies D1–D5 (highest Fe, lowest Mn; Figs. S9 and S10). The Ca and Ca∕Ti profiles generally follow the TIC curve (Figs. 4, S9, S10 and S11b) – less in LT2 due to the high Ti values, supporting a stronger influence of clastic processes – with higher values in unit 5 than in the upper units (Figs. S6b, S10d, e). Br shows a significant second-order relationship with TOC ($R_{\mathrm{bootstraping}}^{\mathrm{2}}=\mathrm{0.42}\pm \mathrm{0.16}$, p<0.001; n=100), although this relationship is stronger for sediment samples with TOC<4 % (Fig. S11f). Br and TOC indicators have the lowest values in banded facies D1 and D2, intermediate values in D3, D4 and D5, and higher values in laminated D6 (Figs. S6c and S10a, b). Facies LT1 and LT2 have intermediate to high Ti, K, Rb, Zr and Si values (Figs. 4, S9 and S10), which are all indicative of minerogenic input to the lake.

4.4 Carbon and nitrogen stable isotopes

The δ¹³C in bulk OM samples ranges between −15 ‰ and −30 ‰, although most values are between −22 ‰ and −28 ‰ (Figs. 4 and S6e), which is in the range of lacustrine phytoplankton and macrophyte OM (Meyers, 2003, 1994). Average values are lower for banded facies (D1–D5) and higher for laminated facies (D6) (Fig. S6e). Values are also higher for LT1 and LT2 sediments, particularly for LT1 ($>-\mathrm{20}$ ‰; Fig. S6e). The δ¹³C versus TOC∕TN ratio boxplot shows higher (less negative) δ¹³C values corresponding to higher TOC∕TN, suggesting that the type of OM (macrophytes versus phytoplankton) is a major control of isotopic values in bulk OM in LdM (Fig. S6e and g).

The δ¹⁵N values range between −1 ‰ and 3 ‰ (Figs. 4 and S6f) with the lowest values in LT1 and LT2, intermediate in laminated facies D6, and highest in banded facies (D1 to D5). Although no isotope data for the modern flora of LdM are available, values obtained from soil, particulate organic matter (POM) in the water, macrophytes and sediments of the nearby El Piojo lake (Fig. 2) range from −1 ‰ to 3 ‰, 14 ‰ to 18 ‰, 2 ‰ to 5 ‰ and −1 ‰ to 2.5 ‰, respectively. Recent sediments in LdM are similar to facies D1 and D2 in cores (Fig. S6), and they show TOC∕TN values between 12 and 14, δ¹³C from −14 ‰ to −20 ‰ and δ¹⁵N from −1 ‰ to 10 ‰, with δ¹⁵N values more positive in offshore than in littoral sediments. The close correlation (Fig. S12) with sedimentary facies and organic composition (algal versus macrophyte) suggests that depositional environments (more littoral versus more distal) play a significant role in N dynamics in LdM.

4.5 Pollen

Pollen spectra are primarily dominated by Poaceae (grasses) and secondly by Ephedra. Poaceae are dominant throughout the record, especially in units 3 and 2, and Ephedra in the upper part of unit 5. The Poaceae ∕ Ephedra ratio shows high values during the first half of unit 5. A major change occurs at 300 cm (Figs. 4, S13), marked by a significant increase in the Poaceae ∕ Ephedra ratio. Although these taxa are not the main component of the vegetation, they are one of the few taxa in the Andean flora with anemophilous pollen dispersion syndrome (high production and dispersal of pollen). Modern pollen rain studies in sites northward and southward of LdM show this altitudinal relationship between Poaceae and Ephedra (Fernández Murillo et al., 2019; Paez et al., 1997). Poaceae pollen is a good indicator of high-altitude vegetation while Ephedra is more abundant in lower vegetation belts; thus an increase in Poaceae suggests vegetation associated with high altitude and more humidity, while Ephedra represents vegetation of lower vegetation belts associated with relatively less humid conditions.

4.6 Age model

We used four approaches to establish the age–depth model of the LdM record and to constrain the recent ¹⁴C reservoir effect: (1) dating the water-dissolved inorganic carbon (DIC) at the mixolimnion (∼20 m); (2) dating modern submerged macrophytes (Table 1); (3) comparing ages from the same stratigraphic interval obtained from macrophytes, terrestrial samples with no reservoir effect (i.e., wood), ²¹⁰Pb and ¹³⁷Cs dates, and tephrochronological markers such as the top ash layer assumed to be from the Quizapú Volcano eruption (1932 CE); and (4) selecting only samples of macrophytes through the sequence for radiocarbon dating (Table 1; Figs. 5, S14). Although the top ash layer has not been geochemically fingerprinted, its estimated age according to ²¹⁰Pb techniques coincides with the Quizapú Volcano eruption (1932 CE), and it can be used as a chronostratigraphic marker in the sequence (Carrevedo et al., 2015). Trace/elemental analysis would be required to definitively establish the provenance of this tephra layer; the available compositional data (XRF, DRX and microscope smear slide observations) and the age mode support our assumption that this layer is the younger Quizapú pinion eruption. According to this age–depth model, the LdM sequence spans the last 13 kyr with a reservoir effect of 4.7±0.6 kyr. The age model is thus more robust for the upper three units (Fig. 5), including ¹⁴C dates from terrestrial macrofossils, ²¹⁰Pb and ¹³⁷Cs dates, the Quizapú ash, and a better constrained ¹⁴C reservoir effect. The 1963 CE ¹³⁷Cs peak centered at 6.5–7 cm (Fig. 5b) fits well with the ²¹⁰Pb chronology and the Quizapú ash horizon at 14–15 cm, adding considerable confidence to the age model for these upper units. Units 1 and 2 span the last 0.65 kyr with an accumulation rate of ∼10 and ∼20 a cm⁻¹, respectively. Unit 3 (195–43 cm) spans between 4.0 and 0.65 ka cal BP with a median rate of ∼23 a cm⁻¹. Unit 4 (LT1; 300–195 cm in site 3) started at 4.0 ka cal BP. Laminated unit 5 (440–300 cm) is characterized by an accumulation rate of ∼70 a cm⁻¹ and includes most of the Holocene (13.0 to 4.0 ka cal BP), whereas unit 6 (LT2) was deposited prior to 13 ka cal BP.

Figure 5 Bayesian chronological model for LdM sequence based on combined analyses of ²¹⁰Pb, ¹³⁷Cs and nine AMS ¹⁴C dates. Bayesian age model of LdM showing the calibrated ¹⁴C dates and the age–depth model (Blaauw and Christen, 2011). (a) Chronological model for the last 14.0 kyr. (b) Detail of the CRS (constant rate of ²¹⁰Pb supply) model, the 1963 depth determined from ¹³⁷Cs peak (red star at 6.5–7 cm) and the Quizapú ash horizon (∼14–15 cm). (c) The accumulation rates (in years per centimeter) as estimated by the MCMC iterations (Blaauw and Christen, 2011) with a median of 10, 20, 23 and ∼70 yr cm⁻¹ for units 1, 2, 3 and 5, respectively. (d) Schematic representation of the volcanic facies ages estimated from the chronological model for LdM. The ash and lapilli layers are color-coded according to macroscopic and microscopic features and compositional data from XRF analysis: pink – ash with higher Si and Rb and lower Fe and Sr; mauve – ash with higher Sr, Ti and Fe and lower Si and Rb; mauve – lapilli layers with MS>200 SI, higher Ca and Sr, and lower Fe and K; pink – lapilli with MS<200 SI, higher K and lower Sr, Ca and Fe. (e) Distribution of the Laguna del Maule eruptive unit ages reported by Andersen et al. (2017): pink (rhyolite), mauve (rhyodacite) and green (andesite) colors represent postglacial eruption units described in Hildreth et al. (2010), Andersen et al. (2017) and Singer et al. (2018) (see Fig. 1). Note in (d, e) the greatest number of volcanic events during the late glacial period–Early Holocene and Middle–Late Holocene transitions.

Figure 6 Seismic-to-core correlation and seismic to lithostratigraphic unit comparison. (a) LdM map showing the location of main seismic lines with line numbers. (b) S–N 3.5 kHz seismic profile crossing LEMA 3 (red line in Fig. 6a). We recognize five different seismic units (S1 to S5), based on the characteristics of seismic reflectors and correlated with lithostratigraphic units in northern areas of LdM. Rectangles mark detailed views shown in Fig. 5c. (c) Seismic-to-core correlation of site 2 along the cross N–S reflection seismic profile LEMA 3. The seismic units correspond to the retrieved sequence, including banded to laminated lake sediments with intercalated mass wasting deposits (transparent to homogeneous seismic facies), volcanic layers (main reflectors) and coarse lacustrine/alluvial facies at the base (see Supplement).

5 Discussion

5.1 Establishing chronological controls to the lacustrine sequence of LdM

Obtaining a robust absolute chronology for the LdM sequence is hampered by the lack of plant macrofossils and a large ¹⁴C reservoir effect (Fig. 5 and Table 1). Large reservoir effects are common in Andean volcanic lakes likely due to degassed magmatic CO₂ (Holdaway et al., 2018; Sulerzhitzky, 1971; Valero-Garcés et al., 1999). Volcanic CO₂ rising from sublacustrine springs is assimilated by submerged and aquatic vegetation, resulting in dates that are apparently too old (Christenson et al., 2015). The DIC-based reservoir effect was similar to that found in living macrophytes sampled at the dam (age ca. 2.4 kyr) but considerably lower than the estimate based on the comparison of samples from the same stratigraphic level: a macrophyte sample (LEM11-3A, 13 cm, ca. 4.8 kyr) and a wood sample (LEM13-3D, 14 cm, ca. 0.85 kyr). This is consistent with the less negative δ¹³C values of the primary production in the lake due to the contributions of mantle CO₂ (less negative than atmospheric CO₂) to the DIC of the lake. This variable range of the reservoir effect underlines the complexity of the carbon cycle in volcanic lakes with likely temporal and spatial variability and also suggests that biological effects on the littoral versus distal environments may be significant. A detailed study of ¹⁴C age variability in modern sediment, aquatic vegetation and organic producers is needed to understand these differences. With the available data and as most dated samples in our sequence were macrophyte remains, we have established a reservoir effect of 4.7 kyr for the sequence resulting from the difference in age among macrophyte and wood samples from the same level (Carrevedo et al., 2015). Analyses of pollen in the LdM sequence show the first appearance of Pinus pollen at 30 cm (ca. 1750–1700 CE), similar to other pollen records in central Chile (Villa-Martínez et al., 2003; Frugone-Álvarez et al., 2017), and the sharp increase in Pinus sp. pollen at 10 cm depth (Table S3) coincides with the promulgations of the forest law in 1931 that had a large impact on new forest plantations (e.g, in El Maule Province they reached 143.450 ha in 1943 CE). We are aware of the simplification of these methods to estimate the reservoir effect for the whole sediment sequence as the complexity of the carbon cycle in volcanic lakes is likely to cause spatial and temporal variability in these values. However, until a more detailed tephrochronology is developed for LdM, this methodology offers the best approach to estimate the age of the sequence. We have considered the uncertainties of the age model in the paleoenvironmental and paleoclimatic implications of the LdM record.

5.2 Depositional dynamic of LdM

According to the genetic classification of volcanic lakes by Christenson et al. (2015), Laguna del Maule corresponds to a caldera lake: the lake was developed in a polygenetic volcanic system (G1), the relationship between volcanic processes and lake formation is strong (R1), the duration of the lake water fill after the eruption is long (T0) and the lake fills a large part of the caldera (L1). The main stages in the depositional evolution of LdM have been characterized based on the sedimentological and geochemical features of the composite sediment sequence for site 3. Seven main depositional phases/events have been identified during the evolution of LdM, corresponding to the emplacement of turbidites LT2 (Phase I) and LT1 (Phase V), deposition in a shallower lake with intense volcanism and high organic productivity (phases II to IV), and a deeper, fluctuating and less productive lake (phases VI and VII).

5.2.1 Phase I: emplacement of turbidite LT2

The ¹⁴C-based chronology of the LdM sedimentary sequence favors a pre-Holocene timing (ca. 13.3 ka cal BP) for the emplacement of LT2. The late glacial–Early Holocene chronological framework based on the LdM lake sequence is older than the one derived from the dating of surface volcanic deposits (Andersen et al., 2017; Singer et al., 2018). Andersen et al. (2017) dated the rle (rhyolite of Loma de Los Espejos) lava at the north end of the lake (Fig. 2) to 19±0.7 ka. This volcanic eruption dammed the lake, raising its level by 200 m and forming a prominent shoreline around the entire basin (Hildreth et al., 2010). The dating of these paleoshoreline outcrops with cosmogenic ³⁶Cl provided minimum age ranges of 9.4±0.4, 8.8±0.6, 7.5±0.3, 6.6±0.6 and 4.2±0.2 kyr ago, suggesting that the lake drained catastrophically in the Early Holocene or later. Using this cosmogenic ³⁶Cl dating, Singer et al. (2018) conclude that the age of this event of the paleoshoreline is ca. 9.4 kyr. This is a much younger age (3500 a) compared to the age provided by our ¹⁴C model (ca. 13 kyr). Further dating in both lake sequence and surface volcanic formations is needed to solve this disparity, since the reservoir effect of ¹⁴C may have been variable over time due to different contributions of volcanic gases to the dissolved gas pool of the lake water.

Sedimentary facies at the base of the recovered sequence do not show evidence for deep hemipelagic deposition. According to our ¹⁴C age model, brecciated volcanic facies in site 1 (unit 7) would be pre-Holocene and could correspond to the volcanic Phase 1 defined by Singer (2014). Alternating hemipelagic and coarse sands with rounded clasts suggest strong fluvial/alluvial transport during a lake phase occurring before the catastrophic pre-Early Holocene drainage of the lake (Singer, 2014; Singer et al., 2018). However, as we see no important depositional changes (i.e., from deeper to shallower facies), hiatuses or erosional surfaces in the lacustrine sequence, the sediments recovered in the LdM sequence could only have been deposited after this catastrophic event.

Emplacement of lacustrine turbidite LT2 could be related to three main triggering processes: (1) intense seismic and volcanic activity (Moernaut et al., 2019) in the basin during the late glacial (¹⁴C age model), (2) intraplate (intraslab or crustal) earthquakes (Van Daele et al., 2019; Wils et al., 2018) related to a postglacial/Holocene uplift (Singer et al., 2018) in combination with tectonic processes of the Troncoso fault, or (3) the drainage event of the lake – if this emplacement was related to the catastrophic drainage of the basin – and in this case it should have occurred later, according to the ³⁶Cl dating of the highest lake shoreline in ca. 9.4 ka cal BP (Singer et al., 2018). The discrepancies between the two age models may only be resolved with additional dating methods (e.g., with individual lipid terrestrial compounds and/or estimating a dynamic reservoir effect through the ¹⁴C DIC record – carbonates) and the development of a tephrochronology for the LdM sequence. Indeed, the occurrence of tephra layers in the LdM sequence provides a unique opportunity for fingerprinting and linking them to the known ⁴⁰Ar∕³⁹Ar-dated and ¹⁴C-dated eruptions in the basin (Andersen et al., 2017; Singer, 2014; Singer et al., 2018).

5.2.2 Phase II: a shallow, low productivity lake

Unit 5 represents deposition in the LdM between two major catastrophic events: LT2 (late glacial to Early Holocene, according to our ¹⁴C age model) and LT1 (mid-Holocene, 4 kyr). After deposition of LT2 turbidite, fine lacustrine deposition resumed in the northern LdM basin with diatomaceous facies D6 (sites 1, 2 and 3; Fig. 3). At site 3, sediments overlying LT2 are only about 30 cm thick (unit 5, 411–438 cm depth) and are composed of facies D6a (less diatomaceous and less organic-rich) and D6b (more organic-rich and with some endogenic calcite). This alternation reflects a transition from diatom-dominated productivity (D6a) to macrophyte-dominated (D6b) productivity, likely caused by lake level fluctuations (Fig. 7). A major volcanic event deposited L5 and several tephra (T18 and T17) (Table 5). Higher Z₁-score values suggest a relatively higher supply of volcanic detritus from the catchment to the lake that will continue during the next phase (Fig. 7). The sharp drop in Fe∕Mn values indicates a rapid change to more oxic conditions at the lake bottom or a large change in water geochemistry as a result of the volcanic activity. Several studies have shown that the Fe∕Mn ratio in the sediment is influenced by the redox conditions at the sediment–water interface at the time of deposition (Frugone-Álvarez et al., 2017; Moreno et al., 2007; Naeher et al., 2013; Wersin et al., 1991). More dominant reducing conditions at the water–sediment interphase would enhance dissolution of Mn (II) but would not greatly affect Fe (II), leading to a higher Fe∕Mn ratio in the sediment (Davison, 1993; Tipping et al., 1981). Furthermore, in a reducing system, the sulfates are reduced to sulfides precipitating FeS, which would increase the Fe∕Mn ratio of the sediment. Differences in the input of Mn and Fe to the sediments are likely to occur depending on weathering intensity and volcanic activity. Higher sediment input from the watershed and more oxic conditions at the bottom of the lake are consistent with lower levels and increased mixing regimes in the lake. Most indicators (TOC, TIC, TOC∕TN, Z₂ score, Br∕Ti) suggest a change in the amount and composition of phytoplankton and terrestrial organic matter (Meyers, 2003, Figs. 7 and S12). At this time, BioSi reached the lowest values of the whole sequence, suggesting lower diatom productivity (Fig. 7). The δ¹³C in bulk OM values is in the same range although slightly more negative than other Andean lakes (e.g., Pueyo et al., 2011; Díaz et al., 2016; Contreras et al., 2018), where sources of heavy carbon attributed to magmatic CO₂ inputs seem to play a significant role in C cycling (Valero-Garcés et al., 1999). The δ¹⁵N values depend on contributions of algae, macrophytes, land plants and the nitrogen source (${\mathrm{NO}}_{\mathrm{3}}^{-}$, ${\mathrm{NH}}_{\mathrm{4}}^{+}$, N₂) (Botrel et al., 2014). δ¹⁵N is also used as an indicator of primary productivity (Meyers, 2003) with less positive values linked to a decrease in phytoplankton productivity (Meyers and Teranes, 2002). Pollen samples indicate sparse vegetation and a relatively high Poaceae ∕ Ephedra ratio, suggesting expansion of the high-altitude belt and retraction of the low-elevation vegetation belt, associated with humid conditions during this period. Therefore, the pollen point to a shift to increased aridity at the base of unit 5 and a return to more humid conditions in upper units 3 to 1 (Figs. 7 and S13).

Table 5 Age estimates for LdM tephra.

Figure 7 Holocene paleoenvironmental evolution in the LdM sequence. Comparison of organic productivity, paleoredox, sediment delivery indicators (Z₁ score) and Poaceae ∕ Ephedra ratio plotted according to the posterior age–depth model (gray). Total organic carbon (TOC, %), opal (BioSi, %), total inorganic carbon (TIC, %), total sulfur (TS, %) and the XRF ratios in logarithm of carbonate (Ca∕Ti) productivity, anoxic conditions (Fe∕Mn), organic (Br∕Ti and Z₂ score) productivity and sediment delivery indicators (Z₁ score). Darker grays indicate more likely calendar ages (age model uncertainty); the white curve shows the single “best” model based on the weighted median age for each proxy; the solid black curve shows the fit for each time series using a smoothing polynomial spline; and orange lines are the average of each time series. The blue, light blue, red and light brown boxes indicate the timing of cold–humid, dry–cold, dry–warm and dry conditions recorded in the LdM, respectively.

5.2.3 Phase III: a carbonate-producing, moderate bioproductivity lake

The dominance of facies D6c indicates a large depositional change in the lake basin. Sediments with a higher content of macrophyte remains, higher TIC and Ca∕Ti values, and the presence of small crystals of endogenic calcite and sulfur-rich minerals (TS, ∼2 % to 7 %; Fig. 7) mark the onset of this phase. The Z₁-score values remained high, indicative of a maintained high detrital input. Although with some peaks, the relatively lower Fe∕Mn ratio suggests dominant oxic conditions at the sediment–water interface (Fig. 7). Sedimentological and geochemical evidence points to an alkaline and oxic environment, which contributed to carbonate formation and sulfur oxidation (Fig. 7). During this phase, the TOC∕TN values show peaks with the highest values in the sequence (∼15); moderate BioSi indicates relatively lower phytoplankton productivity (Fig. 7). Fluctuations in TOC, TOC∕TN, δ¹⁵N, BioSi, Br∕Ti and Z₂ score point to periods of improved environmental conditions for macrophyte growth (more available littoral settings), reflecting a climate shift towards warmer conditions. The close relationship between sedimentary facies and organic composition (phytoplankton versus macrophytes) suggests that the depositional environment (more littoral versus more distal) plays a significant role in C and N dynamics at LdM (Fig. S11). The decrease in TOC and the less abundant calcite occurrences indicate a less productive environment towards the end of this phase. Simultaneously, the decrease in Poaceae ∕ Ephedra values is indicative of increased aridity (Fig. 7).

Carbonate deposition in another Andean volcanic lake – Chungará Lake, 18^∘ S, 4500 m a.s.l. (Fig. 1) – has been extensively studied (Pueyo et al., 2011; Sáez et al., 2007; Moreno et al., 2007) and provides a comparison for this site. LdM is an alkaline volcanic lake (Pecoraino et al., 2015) with pH>8, relatively high ${\mathrm{HCO}}_{\mathrm{3}}^{-}$ and ${\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}$ concentrations, and low Ca²⁺ concentrations. Several factors control carbonate formation in volcanic Andean lakes: (i) variations in salinity due to evaporation, (ii) the input of calcium due to weathering of new volcanic material, (iii) CO₂ photosynthetic depletions related to seasonal phytoplankton blooms, (iv) the development of littoral settings more favorable to charophyte growth and (v) changes in C cycling due to volcanic activity. The source of Ca in LdM is likely the weathering of the andesitic/basaltic and rhyolitic rocks from either andesitic/basaltic lithologies of the older Cola de Zorro Formation or late glacial–Early Holocene rhyolitic eruptions (Hildreth et al., 2010). The increased precipitation of calcite and higher presence of sulfate minerals during this phase could be related to the synergistic effects of increased volcanic activity – demonstrated by the higher number of tephra layers intercalated in the sequence – and favorable environmental conditions conducive towards increased alteration ratios (more glassy and highly porous pyroclasts), which increases alkalinity and releases Ca (and other nutrients) into the lake (Figs. 7 and 8). The increased volcanic activity during Phase III could correspond to the Early Holocene volcanic phase defined by Andersen et al. (2017). Higher temperatures during the Early Holocene could have also been a significant factor for promoting the chemical weathering of the surrounding volcanic rocks, while in Chungará Lake the role of salinity due to the aridity could be a more decisive factor.

Figure 8 Comparison of selected proxies in the LdM sequence with representative local, regional and global paleoclimate records. (a) Dome C (δD) ice core record (Jouzel et al., 2007). (b) Southern Hemisphere summer (W m⁻²) insolation to 35^∘ S (Berger and Loutre, 1991). (c) Iron contents of core GeoB 3313-1 (Lamy et al., 2002). (d, e) Carbonate (Ca∕Ti) and total organic productivity (Br∕Ti) of the LdM sequence; shading indicates age model uncertainty. (f) Sea surface temperatures (SST_alkenone) from the continental slope off mid-latitude Chile (Kim et al., 2002) at ∼33^∘ S. (g–n) Categorization of the moist conditions interpreted by Holocene paleoclimate records in central Chile, where 1 is wet, 0 is similar to present and −1 is dry. (g) Tagua Tagua Lake (Valero-Garcés et al., 1999). (h, i) Aculeo Lake Jenny et al., 2002, 2003. (j) Palo Colorado (Maldonado and Villagrán, 2006). (k) Ñague Swamp Forest (Maldonado and Villagrán, 2006). (l) Quereo Swamp Forest (Villagrán and Varela, 1990; Villa-Martínez and Villagrán, 1997). (m) Negro Francisco Lake (Grosjean et al., 1997). (n, o) Holocene marine record, ENSO-sensitive, from Peruvian shelf (Rein et al., 2005). (p) Bulk Ti content of Cariaco Basin sediment (Haug et al., 2001) (see Fig. 1). ACR – Antarctic Cold Reversal; YD – Younger Dryas. The blue, light blue, red and light brown boxes indicate the timing of cold–humid, dry–cold, dry–warm and dry conditions recorded in the LdM, respectively.

5.2.4 Phase IV: A high carbonate, more productive and macrophyte-dominated lake

Organic-rich, carbonate-bearing facies D6a dominated deposition in the central areas of LdM during this phase. Bioproductivity proxies (TOC, Br∕Ti and Z₂ score) show a trend with elevated values at the beginning of this phase followed by a decreasing up-core trend (Fig. 7). Moderate BioSi values and decreasing TOC∕TN reaching the lowest values in the sequence suggest a period of phytoplankton dominance as organic producers during the beginning of this phase. Dominant oxic conditions continued and together with a decreased detrital input (lower Z₁ score) and lower Poaceae ∕ Ephedra values also indicate water mixing, lower runoff and drier conditions, respectively (Fig. 7). This phase ended with higher detrital input (higher Z₁ score) and the onset of an increasing trend in Fe∕Mn, consistent with an increase in water depth and the establishment of more distal environments at the coring site. The ¹⁴C-based age model suggests this phase would correspond to the mid-Holocene (8–4 ka cal BP). Geochemical proxies underscore the millennium prior to deposition of LT1 (ca. 5–4 ka cal BP) as one with the lowest bioproductivity and the highest clastic input during the Holocene. Higher clastic input could be related to the increased volcanic activity at the end of unit 5.

5.2.5 Phase V: emplacement of turbidite LT1 (ca. 4.0 ka cal BP)

A complex sequence of events at LdM occurred around 4.0 ka cal BP, starting with intense and prolonged volcanic activity and deposition of a thick, multistory volcanic unit including T9–L3–T8, followed by a thick turbidite unit (LT1) and capped by another tephra unit (T7) (Figs. 3 and 4). This volcanic event is unique in the sequence as basal tephra T9 shows a fine-grained, coarsening upward texture and convoluted lamination that could indicate pyroclastic flow processes or transport to the lake by runoff (Fontijn et al., 2016). This 4.0 ka cal BP volcanic event could be linked to activity from the southeast volcanic centers in the LdM basin during Phase 2, which started with the development of the Barranca sequence (6.4 ka cal BP) to the southeast of the lake (Andersen et al., 2017). The thick turbidite LT1 was emplaced in the northern, deeper areas of the basin after deposition of ashfall (T8) and capped by T7 tephra. Depositional mechanisms would have been similar to Early Holocene turbidite LT2. The highest δ¹³C values for organic matter (Fig. 4) could be related to a contribution of volcanic gases to the lake. This deep CO₂ could increase the δ¹³C values of the DIC and the primary production (Valero-Garcés et al., 1999).

5.2.6 Phase VI: a deeper lake during the Late Holocene (4.0–0.7 ka cal BP)

Phase VI started after the emplacement of LT1 with a completely new set of facies (banded D3, D4 and D5) characterized by a sharp decrease in Z₁-score values reaching the lowest values in the sequence and an increase in productivity (TOC, Br∕Ti, BioSi and Z₂ score). The lower organic content in diatom-rich facies D4 compared to macrophyte-rich facies D3 and D5 suggests that the littoral setting is the most productive in LdM, as shown in many other mountain lakes (Michelutti et al., 2015; Quayle et al., 2002; Smol and Douglas, 2007). Fe∕Mn values also showed a sharp increase and maintained higher values until the top of the sequence (Fig. 3). Lower clastic input, more anoxic conditions, higher phytoplankton productivity, higher Poaceae∕Ephedra values, lower δ¹³C values and an increase in facies variability at ca. 4.0 ka cal BP are consistent with more distal/deeper conditions at the coring site, higher/fluctuating lake levels and increased winter precipitation. Modifications in basin morphology and accommodation space after the 4.0 ka cal BP event would have also contributed to this depositional change. Sediment delivery (Z₁ score) increased starting at 1.5 ka cal BP and exhibits centennial-scale oscillations (Fig. 7). A higher frequency of TOC, Br∕Ti, Z₂ score and BioSi values suggest that bioproductivity variability has increased since 4.0 ka cal BP, reaching similar or even higher values than during some mid-Holocene intervals. This increasing variability could be driven by an ENSO-like forcing (Moy et al., 2002). Most proxies (TOC, TOC∕TN, BioSi and Z₂ score) identify two main productivity transitions at ∼2.0 and 0.7 ka cal BP (Fig. 7). The period 4.0–2.0 ka cal BP is characterized by intermediate TOC values and high diatom productivity (BioSi) with some TIC peaks (periods of increased alkalinity); conversely from ∼2.0 to 0.7 ka cal BP, TOC was higher, but BioSi, Ca∕Ti, TOC∕TN, Br∕Ti and Z₂-score values decreased (Fig. 7). Higher Fe∕Mn ratio variability occurred during the first period compared to the second (Fig. 7). Sediment delivery (Z₁ score) increased starting at 1.5 ka cal BP, and exhibits a centennial-scale pattern (Fig. 7).

5.2.7 Phase VII: lower productivity and more anoxic conditions (last 7 centuries)

The depositional evolution of LdM during the last 0.7 kyr has been described in detail in Carrevedo et al. (2015). The onset of this phase is marked by a decrease in bioproductivity and an increase in Fe∕Mn values. After deposition of coarser, more littoral facies D3 during unit 2 (∼0.62 ka cal BP to the late 19th century), finer, less organic-rich facies D1 and D2 were deposited during the last century (unit 1). Bioproductivity indicators (TOC, BioSi, Br∕Ti and Z₂ score) show two century-scale peaks from 1300 to 1400 CE and 1650 to 1850 CE and low values from 1400 to 1650 CE, as well as during most of the late 19th and 20th centuries. Pinus and Rumex pollen at the top of the sequence reflect recent human activity (Carrevedo et al., 2015). The damming of the lake in the mid-20th century increased the average lake level, favoring finer sedimentation and more anoxic conditions at the bottom of the lake. The sedimentological and geochemical proxies are consistent with a decrease in reservoir level during the last decade due to lower water resources in the area.

5.3 Paleoclimate implications

The LdM sequence provides new insights on central Chilean paleoclimatic history, particularly about the nature of the Early Holocene, the timing of the mid-Holocene transition on a regional scale and the nature of the larger climate variability during the Late Holocene. In the following sections we will take into account the chronological uncertainties of the LdM age model when discussing the events in LdM and the paleoclimatic inferences.

5.3.1 A relatively arid Early Holocene

After the emplacement of LT2 (Phase I) and likely after the catastrophic drainage of the LdM basin (Singer et al., 2018), the lake was characterized by the dominance of shallower environments with relatively low productivity (Phase II). Soon afterwards, a large depositional change towards a more carbonate-producing, organic-rich and matter-accumulating lacustrine system occurred with the expansion of environments adequate for macrophyte growth (littoral areas) and a more dominant oxic lake bottom with well-mixed waters. During this phase, bioproductivity indicators display oscillations on a millennial scale embedded within a long-term decreasing trend (Figs. 7 and 8d, e). These cycles suggest relatively lower lake levels with more macrophyte littoral productivity alternating with higher lake levels and increased planktonic productivity (Fig. 7). Pollen ratios (Fig. 7) show high values during the first half of unit 5, between ca. 13 and 10 ka cal BP, suggesting relatively humid conditions. Afterwards, a decreasing trend in the pollen ratio suggests an upward drift of the low vegetation belt associated with drier conditions. These changes have been interpreted as the onset of the Holocene, regionally characterized by the transition from a dry–cold climate during the late glacial to a dry–warm climate during the Early Holocene (Kim et al., 2002). After ca. 9.0 ka cal BP, a slight reversion in trends suggests relatively more humid conditions until ca. 8.0–7.5 ka cal BP. A sharp increase in organic productivity and more abundant macrophytes (TOC∕TN values ∼17, and δ¹⁵N values between 0.5 ‰ and 1 ‰) are indicative of lower lake levels in LdM during the Early Holocene, which could have been caused by a decrease in snow accumulation during winter. High clastic input (higher Z₁ score; Fig. 7) during this phase could be explained by the dominance of shallower, littoral settings in the lake. These lower lake levels, and relatively more arid conditions in the LdM record during the Early Holocene, are in good agreement with other regional records from the continent (Tagua Tagua, about 150 km north of LdM; Valero-Garcés et al., 2005; Jara and Moreno, 2014; Figs. 1 and 8) and the ocean (Lamy et al., 2010, 1999; Muratli et al., 2010a, b). Arid conditions during the Early to mid-Holocene are also shown in pollen records at Quintero and Quereo (33^∘ S) (Villagrán and Varela, 1990) and in Laguna Aculeo (34^∘ S) (Jenny et al., 2002) from >9.5 to 5.7 ka cal BP (Figs. 1 and 8). Drier conditions during the Early to mid-Holocene have also been documented at 40–43^∘ S (Abarzúa et al., 2004; Moreno, 2004; Moreno and León, 2003). Sedimentary records from oceanic, fjord and lake sites from southern, central and northern Chile indicate widespread warming at Southern Hemisphere mid-latitudes between 12.5 and 8.5 ka cal BP (Lamy and Kaiser, 2009; Fletcher and Moreno, 2011, 2012; Kaiser et al., 2008; Kim et al., 2002; Lamy et al., 2010).

The abrupt shift to warmer conditions at the onset of the Holocene (Figs. 7 and 8) could have been caused by synergetic interactions between the southward displacement of the SWW and a strengthening of ENSO, associated with a weakening in the Atlantic Meridional Overturning Circulation (AMOC) that would lead to a southward shift in the Intertropical Convergence Zone (ITCZ) (Haug et al., 2001; Kienast et al., 2006; Mosblech et al., 2012; Rein et al., 2005; Schneider et al., 2014). During the Early Holocene, the SWW would have been in a more summer-like condition compared to a more winter-like pattern during the Late Holocene (Lamy et al., 2010). A positioning of the SPSH at the latitude of LdM (∼36^∘ S) with an initial weakening in the Early Holocene and a subsequent strengthening towards the Late Holocene would have caused SWW-summer conditions between ∼12.0 and 9.0 ka cal BP and decreased moisture in the Andes of central Chile (Figs. 7 and 8).

5.3.2 A complex mid-Holocene transition (8.0–4.0 ka cal BP)

Phase III likely spans the mid-Holocene (8–4 ka cal BP) and represents the period with the highest organic productivity in the lake. The lowest values in the pollen ratio occurred mainly around 7.5 and 5.0 ka cal BP, suggesting the maximum expansion of the lower vegetation belts and drier conditions particularly between ca. 7.0 and 6.0 ka cal BP, the driest period of the entire record. Taking into account all proxies analyzed, productivity started to decrease at ca. 6.0 ka cal BP reaching minimum values around 4.5–4.0 ka cal BP (Figs. 7 and 8). This shift is accompanied by more phytoplanktonic productivity, and it could reflect an increase in water depth between ca. 6.0 and 4.0 ka cal BP as littoral, macrophyte-rich environments were reduced in surface area as shown in other lakes (Nõges, 2009). Hydrological changes beginning in ca. 8–6.5 ka cal BP are also known from other regional records (Fig. 8). Reconstructed precipitation from sedimentological, palynological and diatom data from Laguna Aculeo (∼34^∘ S, 300 m a.s.l.) suggests more arid (150–300 mm a⁻¹) and warmer conditions prior to 7.5 ka cal BP, a progressive increase in effective moisture occurring around 5.7 ka cal BP with the development of a fresher water lake (up to 450 mm a⁻¹), and a marked lake level rise after 3.2 ka cal BP (up to 550 mm a⁻¹) (Jenny et al., 2002). Geochemical and pollen records from the north-central Chilean coastal areas show an open landscape between 8.7 and 5.7 ka cal BP as a reflection of drier conditions, followed by an increase in arboreal pollen during increasingly wetter conditions between 5.7 and 3.0 ka cal BP (Frugone-Álvarez et al., 2017; Maldonado and Villagrán, 2006; Villagrán and Varela, 1990). Low lake levels and high littoral bioproductivity at LdM between ca. 8 and 6 ka cal BP (Fig. 8) are also coincident with decreased winter precipitation in the Andes during the mid-Holocene (Fornace et al., 2014; Moreno et al., 2007; Valero-Garcés et al., 1999, 1996). In addition, a neoglacial advance in the Mendoza Andes valley located at ∼35^∘ S suggests cold conditions between 5.7 and 4.0 ka cal BP (Espizua, 2005), which is coeval with the changes that took place at LdM during this time interval (Fig. 8).

Further north, the Atacama and the tropical Andes records show a complex humidity spatiotemporal structure during the Early to mid-Holocene with clear differences between high-altitude environments (Giralt et al., 2008; Grosjean and Núñez, 1994; Grosjean et al., 2003; Moreno et al., 2007; Sáez et al., 2007; Schittek et al., 2015; Valero-Garcés et al., 1999, 1996) and the mid-altitude areas (Latorre et al., 2003; Betancourt et al., 2000; Quade et al., 2008). Intense and rapid moisture fluctuations during the mid-Holocene associated with modifications in precipitation seasonality have been proposed as an explanation for these discrepancies (Grosjean et al., 2003). South of LdM, at 46^∘ S, increasing wind and precipitation (Daele et al., 2016) occurred between ca. 10 and 5–4 ka cal BP; in southern Patagonia (53^∘ S) the main decrease in precipitation also occurred between 11 and 8 ka cal BP, but this decrease continued until ∼5 ka cal BP (Lamy et al., 2010). The lag in the climatic response from south to north is mainly a result of the expanding/shifting SWW belt toward the north.

Instrumental rainfall records from central Chile (27–36^∘ S) indicate a tendency towards warm–wet (cold–dry) climate during positive (negative) phases of ENSO–PDO that have a strong influence on the snow pack accumulation over the subtropical and mid-latitude Andes (Garreaud and Falvey, 2009; Masiokas et al., 2006; Montecinos and Aceituno, 2003). Variation in the strength of the circulation associated with ENSO–PDO-like teleconnection (e.g., cold fronts that are controlled atmospheres of the SWW and SPSH) could explain a great part of the variability in bioproductivity in LdM: during periods of strengthened ENSO–PDO-like conditions, larger winter snow accumulation could produce higher runoff, increased lake levels, decreased littoral (macrophyte) productivity and increased planktonic productivity. On the contrary, periods of weakening ENSO–PDO-like conditions (or more frequent La Niña-like conditions) would have favored increased organic productivity associated with a lower water level and less intense winters (Masiokas et al., 2006; Nõges, 2009).

5.3.3 A more humid Late Holocene (last 4 kyr)

Increased bioproductivity, decreased clastic input and the onset of dominant anoxic conditions (higher Fe∕Mn) in LdM at ca. 4.0 ka cal BP (Fig. 7) point to higher lake levels and increased phytoplanktonic productivity. The Fe∕Mn ratio shows a change in oxidation state in the sediment–water interface at ∼4.5 ka cal BP to more anoxic conditions that continued throughout the Late Holocene (Fig. 7). Although geomorphic changes in the lake basin due to volcanic/seismic activity might have played a significant role, the consistent presence of high lake levels in other regional records favors a climate influence with higher winter snowfall compared to the Early and mid-Holocene. All bioproductivity proxies show a distinctive transition at 2.0 ka cal BP with a decreasing productivity and a relative increase in clastic input (Fig. 7). Furthermore, the Late Holocene in LdM is characterized by high pollen ratio values, suggesting the expansion of the high-altitude belt and retraction of the low-elevation vegetation belt, associated with humid conditions during this period.

Overall, moisture availability in central Chile during the last ca. 4.0 kyr seems to have increased, paralleling summer insolation and seasonality, which reached their maximum in the Late Holocene (Wanner et al., 2008). The Aculeo record shows a progressive increase in effective moisture after 5.7 ka cal BP and the establishment of modern humid conditions around 3.2 ka cal BP (Jenny et al., 2002; Villa-Martínez et al., 2003). Numerous intercalated clastic layers after 3.2 ka cal BP reflect floods during rainy winters associated with more frequent/intense El Niño-like activity (Jenny et al., 2002). Furthermore, paleoceanographic proxy data from a marine core off central Chile show a significant increase in paleoproductivity and decreasing SSTs at 33 and 41^∘ S (Kim et al., 2002; Lamy et al., 1999), which is also interpreted as increased ENSO activity beginning in ca. 5.0 ka cal BP (Rodbell et al., 1999). Higher lake levels and increased moisture availability during the Late Holocene have also been documented in the Altiplano (Laguna Miscanti, 23^∘ S; Valero-Garcés et al., 1999, 1996; Grosjean et al., 2001) and in northwest Patagonia (Lago Pichilafquén, ∼40^∘ S; Jara and Moreno, 2014). The LdM record supports increasing moisture during the Late Holocene from tropical to mid-latitudes in South America. Both greater fluctuations in water levels and stronger seasonality in temperature and precipitation would have prolonged seasonal cycles and modified the hydrological dynamics of the LdM basin. The invoked climate mechanisms include a strengthening of ENSO–PDO-like activity and a northward shift of the ITCZ (Bird et al., 2011; Carré et al., 2014; Haug et al., 2001; Moy et al., 2002; Rein et al., 2005; Stansell et al., 2013; Zhang et al., 2014) (Fig. 8).

The LdM productivity record shows two main changes during the last millennium: (i) around 1300 CE, when a drop in organic bioproductivity and the onset of more frequent anoxic conditions occur, and (ii) around the late 19th century. The first period matches with changes in the southern Pacific, leading to a cooler and wetter climate with increased storminess due to ENSO intensification, which largely affected the Polynesian societies (Rull et al., 2015). The second, and more significant, fall in productivity is concomitant with the rapid increase in Aulacoseira and the decrease in planktonic Discostella stelligera and marks the end of the Little Ice Age (LIA; 600–150 a cal BP) (Carrevedo et al., 2015). Pollen data indicate relatively colder and moister conditions during the LIA (Carrevedo et al., 2015). The Medieval Climate Anomaly is characterized by a higher bioproductivity in comparison to the LIA (Figs. 7 and 8).

5.4 Holocene volcanic history based on LdM tephra deposits

The LdM sequence offers a continuous record of volcanic activity in the basin (Fig. 5d), which can be compared to the sequence described previously by Andersen et al. (2017) on surface deposits (Fig. 5e). Nevertheless, the assignment of these events found in the lake sequence (Fig. 3) to dated volcanic eruption (Singer, 2014; Singer et al., 2018) is not straightforward. Detailed mapping and dating of surface volcanic formations (Fig. 2) have led to the identification of two main recent phases of silicic volcanism in LdM: early postglacial (between 26 and 19 kyr) and Holocene (Andersen et al., 2017; Singer, 2014; Singer et al., 2018). The occurrence of tephras T23 to T19 is consistent with regional volcanic activity prior to the Holocene (Fig. 3). The presence of coarser tephras (Lapilli L6, site 1) suggests activity from closer volcanic centers. Although a genetic relationship could not be demonstrated, the relationship with LT2 suggests that the emplacement of L6 occurred at the same time as major volcanic/seismic events in the LdM basin (Fig. 3). Intense shaking of slope sediments and instability caused by hydrological changes or sediment load related to volcanic activity could have triggered mass wasting processes resulting in the LT2 emplacement in the basin (Fig. 3). This is consistent with the lower δ¹³C and TOC∕TN ratio values which are indicative of the greater content of C₃ land plants in the sediment (Fig. S12). Central Chile is also well-known for its megathrust earthquakes, and many Andean lakes contain turbidites triggered by earthquake shaking (Van Daele et al., 2015). Turbidites like LT1 and LT2 could be triggered by intense volcanic activity and/or large megathrust earthquakes that caused the deposit of a thick, coarse-grained tephra in the lake, lowering the slope-failure threshold (Wiemer et al., 2015).

Volcanic activity in the LdM volcanic field during the Early Holocene was intense (Fig. 3 and Table 5). After the large volcanic event responsible for L5, four tephra layers were deposited (T17 to T14). More silica-rich composition of tephra (T18 till T9) and lower MS values of L4 and L3 suggest a dominant rhyolitic chemistry composition during this phase, and it is also consistent with the end of the Late Pleistocene to Holocene volcanic phase (Hildreth et al., 2010; Singer et al., 2000; Singer, 2014; Andersen et al., 2017). Some tephra from the Early Holocene (T14, T15 and T16; see below) also have a distinctive composition with low XRF values of Fe and Sr (Fig. S7). Interestingly, this was a period with high organic content in the sediments (more macrophyte remains in facies D6c), the highest TS and the common presence of endogenic calcite, all suggesting some synergetic interactions between volcanic activity and carbonate formation in the lake (Fig. 4). The second half of this volcanic phase was relatively less active and ended with the deposition of L4 (Fig. 5). Although the dating of this part of the sequence is still uncertain, the onset of calcite formation in the lake could be synchronous to the beginning of the Holocene, and so this phase could encompass the Early Holocene (up to 8 ka cal BP).

Lapilli L2 and L1 and ash T6 to T1 layers deposited during the last 4 kyr have petrographic and chemical characteristics, pointing to a silica-poor composition possibly associated with basaltic-andesitic magma (Fig. S7). Eruptive volcanism in LdM during the last 4 kyr issued from eruptive centers around the areas of Cari Launa (rsl: 3.3±1.2 kyr; rcl: <3.3 kyr), Divisoria (rcd: 2.1±1.3 kyr), Laguna Sin Puerto (rdsp: <3.5 kyr), rhyodacite of the northwest coulee (rdcn; 3.5±2.3 kyr) and Colada Las Nieblas (rln: Late Holocene) (Figs. 2 and 5). Eruptions during the Holocene phase are dominantly rhyolitic eruptions, although andesitic, rhyodacitic and mixed eruptions also occurred (Figs. 2 and 5). However, the Late Holocene rhyodacite eruptions of Laguna Sin Puerto and northwest coulee could be linked to the deposition of L2 and L1 lapilli.

6 Conclusions

The composite LdM sediment sequence includes distal lacustrine, volcanic and massive wasting deposits. Six lithostraphic units have been defined in the northern area of the basin and correlated with five seismic units. Lacustrine Turbidite LT2 (unit 6) is composed of massive black silt facies, whereas LT1 sediments (unit 4) are browner, coarser material with abundant macrophyte remains and also millimeter-size pumice clasts. Biogeochemical differences between LT1 and LT2 could imply different depositional processes (seismic and volcanic) and/or provenance. Although some discrepancies between our age model and the dating of volcanic episodes concerning the timing of some major hydrological events (blockage of the outlet by lava flows, highest lake level indicated by shorelines around the lake and drainage of the lake) remain, the composite LdM sequence spans the Holocene after the catastrophic drainage of the lake basin likely due to upstream erosion of the Maule River. Volcanic facies occur as lapilli (6 layers) and ash (23 layers). Their compositional features suggest a Late Holocene transition towards a more silica-rich magma composition. Despite the chronologic uncertainties, the LdM record indicates lower lake levels during the Early Holocene with millennial-scale bioproductivity changes consistent with lower summer insolation and increased aridity. Higher bioproductivity occurred during the mid-Holocene (from ca. 8.0 to 6.0 ka cal BP), synchronous with the phase of aridity described for the tropical and temperate latitudes of South America. During the Middle to Late Holocene, the LdM record indicates relatively higher lake levels, consistent with increased moisture after 4.0–3.0 ka cal BP, caused by the inception of the current ENSO–PDO-like dynamics in central Chile. The Little Ice Age shows a two-phase structure with cold/wet intervals between 1300 and 1450 CE and 1600 and 1850 CE, interrupted by a warmer climate between 1450 and 1600 CE. The LdM record also suggests that millennial-scale Holocene climate and water availability in central Chile were largely ruled by variations in the summer insolation. Complex interrelations between solar irradiance and dynamic changes in regional patterns of internal climate variability, such as the ENSO–PDO, SWW and the SPA, however, seem to exert major control on centennial to decadal scales.