High-resolution reconstruction of infiltration in the Southern Cook Islands based on trace elements in speleothems

Abstract

This study utilizes speleothem trace elements as climate proxies to reconstruct hydroclimate variability over approximately 350 years in the Southern Cook Islands. Stalagmites Pu17 and Pu4 from Pouatea cave were analyzed using high-resolution LA-ICP-MS for trace elements (Mg, Na, Sr, P, U, Y). By monitoring cave dripwater and conducting regression analysis, we found that Mg, Sr, and Na in Pouatea dripwater mostly originated from marine aerosols, while Sr and Ba were primarily from bedrock, with additional Ba coming from marine aerosols and weathered oceanic basalt leaching. Mg was identified as the most reliable element for hydroclimate reconstruction due to its predominantly marine aerosol origin. Infiltration, via dilution of marine aerosols and bedrock inputs, was identified as the main driver of trace element variations in Pouatea at a seasonal scale. Transfer functions were established between each trace element and effective infiltration was calculated, with Mg showing the strongest correlation. The reconstructed infiltration data were compared with climate indices, showing an overarching role of the SPCZ and ENSO in controlling rainfall in the South Pacific. This research demonstrates the potential of speleothem trace elements for paleohydroclimate reconstructions, improving understanding of rainfall variability in the climatically vulnerable South Pacific Islands over the past millennia.

INTRODUCTION

The people of the South Pacific islands are highly vulnerable to the effects of both climate variability and climate change. The region is affected by climate extremes, leading to floods and droughts, with acute effects on agriculture (Held and Soden, 2006; Xie et al., 2010; Barnett, 2011; Cai et al., 2015). Climate variability in the South Pacific is influenced by seasonal shifts of the South Pacific Convergence Zone (SPCZ): the most significant rain belt in the Southern Hemisphere and a key driver of rainfall variability in the region (Higgins et al., 2020). The SPCZ affects Pacific island communities by modulating hydroclimate extremes (Brown et al., 2020), which collectively have affected more than 3.4 million people in the Pacific region and caused over 10,000 reported fatalities since 1950 CE (World Bank, 2017).

Interannual shifts of the SPCZ to a more zonal position and/or its contraction are associated with changes in El Niño Southern Oscillation (ENSO), which then influences the strength of the subtropical anticyclone of the East Pacific and the transport of moisture towards the equator (Trenberth, 1976; van der Wiel et al., 2016). During very strong El Niño events, the SPCZ shifts towards the equator, merging with the Intertropical Convergence Zone (ITCZ) and having a near-zonal orientation (Vincent et al., 2011). On longer timescales, the Interdecadal Pacific Oscillation (IPO) influences both the expansion/contraction of the SPCZ and the sequencing of ENSO events (Salinger et al., 2001; Kumar et al., 2006; Verdon and Franks, 2006; Weir et al., 2021). An analysis of tropical South Pacific rainfall over the period 1961–2000, based on daily rainfall records from 22 Pacific stations, including the Cook Islands, demonstrated that the SPCZ position influences both mean monthly precipitation and daily rainfall extremes (Griffiths et al., 2003). Although much progress has been made towards understanding the SPCZ links between rainfall trends, rain intensity, episodic drought, and deluge, these assessments are increasingly questionable against a backdrop of contemporary warming. What is needed, therefore, are long, quantitative records that capture greater variability in SPCZ dynamics (Magee and Verdon-Kidd, 2019). The limited availability of climate observations in the South Pacific is sufficient justification for obtaining paleoclimate reconstructions from natural archives, such as speleothems, for understanding the past variability of the SPCZ (Brown et al., 2020).

Speleothems provide an increasingly popular archive for terrestrial paleoclimate reconstruction on various time scales (Fairchild and Baker, 2012) and in different climate settings (Bar-Matthews et al., 1999; Li et al., 2005; Wang et al., 2019; Vansteenberge et al., 2020). Speleothems contain both physical and chemical multi-proxy data on surface processes that can reveal causes of climate change on timescales from decades to millennia to millions of years before present (Woodhead et al., 2010). Physical speleothem proxies are related to distinct layering in calcium carbonate mineralogy and/or fabrics, which usually appear as seasonal alternations between dark, compact laminae and white, porous laminae (Railsback et al., 1994; Genty and Quinif, 1996). The origin of visible layering is attributed to changes in the geometric arrangement of crystals, porosity, and/or variability in the concentration of particulates that result from changes in hydrology or colloidal inorganic/organic matter input into the cave (Borsato et al., 2021).

Chemical speleothem proxies include ratios of stable isotopes of oxygen (¹⁸O/¹⁶O, expressed as δ¹⁸O) and carbon (¹³C/¹²C, expressed as δ¹³C), and the concentration of trace elements. Speleothem carbon isotope variations (δ¹³C) reflect soil microbial activity, root respiration, vegetation type above the cave, and rainfall amount (Genty et al., 2001; Johnson et al., 2006; Frisia et al., 2011; Fohlmeister et al., 2020). Speleothem δ¹⁸O signals record variations in cave temperature and local rainfall (e.g., McDermott, 2004; Fairchild and Baker, 2012), and are the most widely applied paleoclimate proxy with which to reconstruct the past (e.g., Wang et al., 2001; Burns et al., 2002; Bar-Matthews et al., 2003; Frisia et al., 2005; Cheng et al., 2006; Mattey et al., 2008; Pausata et al., 2011; Chakraborty et al., 2016).

In addition to commonly used stable isotopes, concentrations of trace elements (e.g., Mg, Sr, Ba, Na, U, P) increasingly are being utilized as tracers of past climate variability (Fairchild et al., 2000; Treble et al., 2003; Borsato et al., 2007; Fairchild and Treble, 2009; Griffiths et al., 2010; Bao et al., 2023; Sliwinski et al., 2023). Ratios of Mg/Ca, and Sr/Ca are often linked to hydroclimate processes via, for example, interpreting their co-variation in terms of rock-water interaction (RWI), prior calcite precipitation (PCP), and drip rate variability (McDonald et al., 2004; Cruz et al., 2007; Sinclair et al., 2012; Treble et al., 2015). Via PCP, Ca ions precipitate as calcite (CaCO₃), resulting in enrichment in Mg and Sr in the remaining solution. Thus, low infiltration conditions (dry periods) result in PCP marked by coupled increases in Mg/Ca and Sr/Ca (Fairchild and Treble, 2009; Belli et al., 2017). Mg and Sr can then be reliably interpreted as proxies of rainfall, particularly when they are employed together with other proxies such as δ¹³C and δ¹⁸O. A positive covariation between δ¹³C and Mg stands out as indicative of a long water residence time within the aquifer during dry periods (Fairchild et al., 2000; Sinclair, 2011; Sinclair et al., 2012). The other hydrologically sensitive elements involve Na, U, and P, which have scarcely been exploited, possibly due to the complexity of their source, transmission through the karst aquifer, and subsequent incorporation in the speleothem. Examples of using these elements as proxies to infer paleohydrology include the use of U as a proxy for water excess (Kaufman et al., 1998; Hellstrom and McCulloch, 2000) and P as a proxy for hydrological and ecosystem variability (Borsato et al., 2007; Lewis et al., 2011; Johnston et al., 2021). Nevertheless, trace element proxies have, so far, mainly been used for qualitative reconstructions of past climate variability and past precipitation, and examples of quantitative reconstructions are limited (e.g., Warken et al., 2018). The dearth of quantitative reconstructions using trace elements largely owes to the complexity of establishing a transfer function between a climate parameter at Earth's surface and speleothem proxy signals. For quantitative reconstructions, a rigorous understanding of the processes involved in proxy signal transfer from the surface into the cave as well as the in-cave processes is needed (Baker and Bradley, 2010; Casteel and Banner, 2015; Warken et al., 2018). Additionally, most of the available reconstructions belong to climate settings where there is a well-emphasized climate seasonality, such as temperate climates (e.g., Baldini et al., 2002; Treble et al., 2003) and monsoon-dominated climates (e.g., Wang et al., 2001). This highlights the necessity of providing an accurate paleoclimate reconstruction for less studied settings such as the tropical South Pacific.

The caves of Atiu Island, Southern Cook Islands, provide excellent material for paleoclimate reconstruction in an area of the South Pacific where the SPCZ dominates. Recent research has shown that stalagmites from Atiu can yield accurate hydroclimate reconstructions from an area that is very sensitive to SPCZ shifts (Faraji et al., 2021). The great advantage of the Atiuan stalagmites is that they have been precisely dated by coupling radiometric ages with annual lamina counting, thus overcoming the common issue of dating difficulties that affect many tropical Pacific speleothems (see Faraji et al., 2021, 2023). Further, caves in Atiu are relatively shallow, with a rock burden of 2–10 m. This ensures rapid transmission of surface-climate conditions into the cave, further enhanced by a limited and patchy soil cover. Except for some small pockets of colluvial red clay, which fill in the joint-controlled karst (Stoddart et al., 1990), the surface above the cave is bare karst. Therefore, speleothems forming in these caves have the potential to encode a climate signal that is rapidly transmitted and, possibly, insubstantially altered.

For this study, we retrieved two stalagmites from Pouatea cave located on the southwestern side of Atiu to reconstruct the variability of effective infiltration (rainfall minus potential evapotranspiration) in the past 350 years by using a cluster of hydrologically sensitive trace elements (Mg, Sr, Na, U, P, Y). This study is one of the first to use speleothem trace elements to achieve high-resolution paleoclimate reconstructions from the South Pacific. We first describe the methods and then present the results and discuss the relationship between the reconstructed infiltration and calculated infiltration. The reliability of our reconstruction is then evaluated against the instrumental records of SPCZ and ENSO.

SAMPLE AND SITE DESCRIPTION

Southern Cook Islands and Atiu Island

The Cook Islands (CI) is an archipelago divided into Northern and Southern groups, situated 9–23°S and 156–167°W in the South Pacific Ocean (Fig. 1b), comprising 15 atolls and coral-rimmed volcanic islands with a total land area of ~238 km² (Merlin, 1985). In the Southern Cook Islands (SCI), easterly trade winds blow for over half the year, mostly oriented 80–140°, depending on the season. During spring, 65% of all winds constantly blow from the easterly quarter, whereas during summer, the SCI are affected by occasional tropical cyclones and extratropical fronts (Thompson, 1986).

Figure 1. (a) South Pacific Convergent Zone (SPCZ) and trade winds in the South Pacific Ocean (redrawn after Australian Bureau of Meteorology and CSIRO, 2011a). (b) The Northern and Southern groups of the Cook Islands (red inset in a; redrawn after Australian Bureau of Meteorology and CSIRO, 2011b). (c) Geological map of Atiu Island (redrawn after Stoddart et al., 1990), with an arbitrary section of raised carbonate rim (A–B section); also showing the location of the Pouatea cave entrance. (d) The surveyed map of Pouatea cave showing the cave entrance, skylights, and locations from which stalagmites Pu4 and Pu17 were retrieved.

Hydroclimate in the SCI is influenced by shifts of the SPCZ (Widlansky et al., 2011; Lorrey et al., 2012; Harvey et al., 2019), which is the largest rain belt in the Southern Hemisphere (Fig. 1a), ENSO (Kruk et al., 2015; Lorrey and Fauchereau, 2018; Weir et al., 2021), IPO, and tropical cyclones (TCs) (Diamond et al., 2012, 2013; Magee et al., 2020). To date, research on the SPCZ indicates that north–south shifts in mean SPCZ position are the key driver of precipitation variability in the SCI, which is positioned at the southern margin of this zone (Rasbury and Aharon, 2006; Partin et al., 2013; Maupin et al., 2014). Shifts in the position of SCPZ appear to be related to both oceanic and atmospheric processes (Widlansky et al., 2011; Brown et al., 2020). The formation of SPCZ is explained through the “graveyard” hypothesis, which suggests that small variations in sea surface temperature (SST) generate upper tropospheric wind fields, which change the wave accumulation patterns. The SPCZ is where these synoptic disturbances tend to accumulate (Widlansky et al., 2011). The SPCZ is variable from month to month (Thompson, 1986; Lorrey et al., 2012), but is an ever-present phenomenon in the South Pacific, thus wielding significant influence on the island communities. The southwest movement of the SPCZ over the SCI results in unsettled weather and heavy rain from November to April, when over 60% of annual precipitation occurs. During the dry season, from May to October, the SPCZ tends to lie farther to the north of the group, and during this time the Southern Cook Islands are affected predominantly by the dry southeasterly trade winds (Thompson, 1986; Stoddart et al., 1990; Barstow and Haug, 1994).

Atiu Island is the third largest of the SCI archipelago (Fig. 1c). The island receives an annual rainfall of 1950 ± 350 mm/year (1986–2018), with 68% occurring during the wet season (November through April). Despite this clear seasonality, rainfall amounts are still considerable during the drier season (May through October). The temperature in Atiu is relatively constant year-round, with a mean of 24.2 ± 1.6°C (1993–2018). The infiltration (water excess = the amount of rainfall minus the amount lost by possible evapotranspiration) was calculated using the Thornthwhaite formula (Thornthwaite, 1948). Monthly rainfall and temperature measurements for Atiu from 2014 through 2018 revealed water excess occurring all year except in September and October. Water excess demonstrates a clear seasonality, with most of the excess occurring from December to May (Fig. 2).

Figure 2. Monthly climate in Atiu. The theoretical potential evapotranspiration calculated from the Thornthwaite equation (Thornthwaite, 1948) based on the monthly temperature recorded above the Pouatea cave from 2014 to 2018.

Geologically, Atiu is a highly eroded volcanic island surrounded by a makatea: a rim of elevated, heavily weathered Cenozoic reef limestone (Fig. 1c) (Stoddart et al., 1990). The volcanic rocks consist of a basalt-basanite suite, ranging from alkali basalt with abundant olivine to olivine-poor, aphyric basalt. When altered, the most abundant secondary minerals are chlorite, clay, calcite, and Fe-oxides. The most common oxides analyzed in the basaltic rock of Atiu are silica, magnesium oxide, calcium oxide, iron oxides, and titanium dioxide (Wood, 1978).

The highly eroded volcanic core reaches a maximum height of ~71 m asl and is mostly covered with limonitic-nodule-bearing red clay. An age of 8–10 Ma is attributed to the volcanic rocks, indicating two phases of eruption (Turner and Jarrard, 1982). The basalt weathering and erosion products may have been transported into the carbonate factory that rimmed the volcanic island. The uplifted reefal limestone that once rimmed the volcanic island has been dated, based on foraminifera content, to Plio-Pleistocene age (Marshall, 1930). Atiu has been thus classified as a composite island with elevated reef limestone (Mylroie and Vacher, 1999). In these composite islands, the volcanic rocks form a catchment surface from which allogenic water is directed to the raised reef carbonate, penetrating deeply into the vadose zone, and thereby maintaining its ability to dissolve the carbonate rock (Mylroie and Vacher, 1999). This process has led to the formation of a maze of cave passages. It is also possible that in such young carbonates, the water flow operates to both dissolve the rock and diagenetically alter it. That is, simultaneously widen voids by solution and transform minerals such as aragonite into calcite, leading to early cementation of voids (Mylroie and Vacher, 1999). Figure 1c shows the geological map of Atiu and a general profile of the raised Cenozoic carbonate rim (Seymour, 1965; Stoddart et al., 1990).

The soil on the makatea is concentrated in dissolution pits and trenches. Except for some small pockets of red, colluvial clay soil filling the bottoms of joint-controlled karst (Stoddart et al., 1990), the surface above the cave is bare karst. Tree roots find their way underground through dissolution cracks and voids, searching for water in the absence of soil. What soil there is consists of Carbonatic Mollisol type and has thin, dark mollic epipedons (soil A horizon) and a very shallow contact with the underlying limestone (Bruce, 1983). Therefore, despite a thick vegetation cover, soil does not cover the whole limestone surface, and the organic layer observed above the caves is mainly of leaf litter. The indigenous forest cover above the caves consists of native Elaeocarpus tonganus and Hernandia moerenhoutiana. Alien species, most typically coconut palms, mark human influence (Holland and Olson, 1989). However, many weedy plants have become naturalized in the central volcanic interiors of Atiu. These ‘aliens’ generally have not been able to spread into the makatea.

Pouatea cave and speleothems

Pouatea cave (20°01'12"S, 158°07'10"W) is a network cave system with several intersecting passages and side galleries and a total surveyed length of 1200 m (Fig. 1d). The cave's main entrance is a vertical shaft with a drop of about 4 m that opens 23 m above sea level and 525 m from the shoreline. There are five other entrances (skylights) with diameters of 3–10 m that were formed by cave roof collapse. The rock burden above the galleries usually is 4–8 m, and is still characterized by high primary porosity, likely due to the relatively young age of the reef and incomplete diagenesis. Such porous epikarst ensures rapid transmission of surface climate conditions into the cave, further enhanced by a limited and patchy soil cover.

Stalagmite Pu17 was collected from Pouatea cave in March 2019 from a gallery leading to the cave's southern end. Pu17 was fed by a relatively slow and constant drip (1 drop every 15 minutes) and was actively growing at a depth of ~7 m beneath the surface when removed from the cave (Faraji et al., 2021). Pu17 is a 53-mm long, candle-shaped stalagmite (Miorandi et al., 2010) with a honey-colored appearance characterized by alternating brown-orange compact translucent layers (CL) and whitish-yellowish porous layers (WPL).

Stalagmite Pu4 was retrieved in October 2014 at a depth of ~7 m beneath the surface towards the cave's southern end (Fig. 1). Pu4 was actively growing at the time of collection and fed by a relatively fast dripping straw with 1 drop every 18 seconds. From October 2014 to September 2018, a thick layer of calcite formed on the drip logger placed over the stump of the removed Pu4, confirming that the stalagmite was active at the time of collection. Pu4 has an incipient candle-shaped morphology with a length of 25 mm and is characterized by a whitish translucent appearance consisting of compact to open columnar fabric (Frisia et al., 2000).

The most common instrument to establish the age of speleothem layers (i.e., ion-counting multi-collector Inductively Coupled Plasma-Mass Spectrometry [ICP-MS]; see Hellstrom, 2003, and references therein) was applied on both Pu17 (Faraji et al., 2021) and Pu4 (Supplementary Information: Table S1, Fig. S1). U–Th dating was carried out on 18 powdered samples (~200 mg) drilled from the stalagmite growth axis, 10 from Pu17 and eight from Pu4. The samples were dissolved in nitric acid, spiked using a ²²⁹Th-²³³U-²³⁶U tracer, and U and Th were eluted in selective ion exchange resin following the procedure in Hellstrom (2003) and Drysdale et al. (2012). The measurements were performed using a Nu Instruments Plasma MC-ICP-MS at the University of Melbourne. The U content was 20–100 ppb, a range observed in active and modern stalagmites from Atiu (Faraji et al., 2023; unpublished materials). An age-depth model was determined from the remaining eight analyses using the finite positive growth rate model of Hendy et al. (2012) and the assumption of active growth at the time of collection.

Tropical island stalagmites, however, can be affected by large dating uncertainties related to their elevated and variable initial ²³⁰Th/²³²Th ratio, which is a commonly reported but as yet poorly understood issue (Beck et al., 2001; Partin et al., 2007; Hua et al., 2012). Annually laminated stalagmites, in this case, can provide good chronological constraints (Baker et al., 2021). This includes acquiring precise relative age by using properties of annual variations in visible growth laminae (Frisia et al., 2003; Tan et al., 2006; Baker et al., 2015), or in the cyclicity of their geochemical properties, such as trace element concentration (Treble et al., 2003; Johnson et al., 2006; Borsato et al., 2007; Orland et al., 2014; Jamieson et al., 2015; Nagra et al., 2017; Ban et al., 2018; Wang et al., 2019), or C and O isotope ratios (Treble et al., 2005a; Mattey et al., 2008). Following the methodology outlined in Faraji et al. (2021), accurate age models were acquired for Pu17 and Pu4 via integrating, in a multivariate analysis, high resolution (6 μm) variations in trace elements that were analyzed by LA-ICP-MS, with optically visible growth bands and two-dimensional Sr-concentration laminae identified through synchrotron-radiation-based micro XRF mapping. These three independent analytical techniques produced more than 30 peak counting series for each stalagmite, 20 of which were averaged and integrated into a single master chronology. The Pu17 chronology was further refined by using the ¹⁴C bomb-pulse peak method (Faraji et al., 2023). The top 12.5 mm of Pu17 was dated using the ¹⁴C bomb-pulse soil continuum method from Markowska et al. (2019) because the very recent laminations are not evident, and fabrics and growth patterns are complicated (Frisia et al., 2022), likely due to the simultaneous dissolution and precipitation of calcium carbonate. The radiocarbon chronology and available laminae-counting age model were then spliced to achieve a single master chronology for Pu17. These combined approaches resulted in age models with a mean 2% uncertainty for Pu17 and 4% uncertainty for Pu4, considerably improving upon the ~50% uncertainty in the U–Th ages (Faraji et al., 2021) (Supplementary Information, Fig. S1).

The Covid-19 outbreak and the remoteness of Atiu allowed only four visits for the purpose of cave monitoring. We monitored Pouatea cave by conducting in situ measurements of pH, electrical conductivity (EC), cave pCO₂, temperature, and humidity. This was accompanied by collection of dripwater samples for analysis. Data loggers were also placed under several dripping straws in the cave's main chamber from 2017 to 2019 to measure temperature, humidity, and drip rate, taking a measurement every hour. The recorded measurements returned a mean air temperature of 23.6 ± 0.5°C for the passage where Pu17 was found, with minimum values in September–October and maximum values in March–April. Relative humidity was near saturation (98–100%) throughout the year, while cave air pCO₂ varied from 500 ppm to 650 ppm in the relatively dry and wet seasons, respectively. Drip rate monitoring indicated a relatively dry period between May and October, which was characterized by a gradual decrease in the dripping rate. The pH of the dripwater in the cave is ~8, and the saturation index for calcite (SIcc) varies from 0.9 to 1. Detailed information about the Pouatea cave monitoring can be found in Faraji et al. (2021).

DATA AND METHODS

Climate data

For this study, we used a range of observed and reconstructed climate indices that represent the climate modes (ENSO and SPCZ) and instrumental precipitation records (Fig. 3). These include local records of monthly precipitation available for four SCI stations (islands): Mangaia (Jan 1914–May 1996), Aitutaki (Jan 1930–Jun 1996), Rarotonga (Jun 1948–Feb 2016), and Atiu (Oct 1980–Mar 2019). We used the instrumental and reconstructed SPCZI from Higgins et al. (2020), which is available at the NOAA Paleoclimate Archive (https://www.ncei.noaa.gov/products/paleoclimatology) (accessed January 2022). SPCZI is a measure of the mean sea level pressure (MSLP) difference between Apia, Samoa (13°48′S, 171°48′W), and Suva, Fiji (18°9′S, 178°26′E), normalized to the 1933–1992 period (Salinger et al., 2014; Higgins et al., 2020), which captures monthly east–west variations in the position of the SPCZ (Higgins et al., 2020). The instrumental record of the Southern Oscillation Index (SOI—representing ENSO) is a normalized measure of the air pressure difference between Tahiti (17°40′S, 149°25′W) and Darwin (12°26′S, 130°50′E) and was acquired from the Australian Bureau of Meteorology (http://www.bom.gov.au/climate/enso/soi/) (accessed January 2022). El Niño events are characterized by a negative SOI (< −5), and La Niña events are characterized by a positive SOI (> +5). The Niño 3.4 Sea Surface Temperature (SST) is another widely used index, which is the anomaly in SSTs in the Niño 3.4 region of the central Pacific (5°N–5°S and 170–120°W). Niño 3.4 SST (based on the Hadley Centre Global Sea Ice and Sea Surface Temperature, HadISST v1.1) was downloaded from the NOAA Paleoclimate Archive, along with the reconstruction of ENSO variability (Li et al., 2011). The instrumental record of IPO was acquired from the Met Office Hadley Centre for Science and Services (https://data.mfe.govt.nz) (accessed January 2022). The two reconstructed records of IPO were accessed from the Australian Antarctic Data Centre (https://data.aad.gov.au) (accessed April 2022) (Vance et al., 2022) and the NOAA Paleoclimate Archive (Linsley et al., 2008). IPO is defined as the difference between the SST anomaly (SSTA) averaged over the central equatorial Pacific and the average of the SSTA in the Northwest and Southwest Pacific (Power et al., 1999). IPO is a large-scale and long-period oscillation that affects a much wider range of latitude than ENSO.

Figure 3. The instrumental and reconstructed climate data used in this study. ENSO = El Niño Southern Oscillation; IPO = Interdecadal Pacific Oscillation (IPO[1] data from Vance et al., 2022; IPO[2] data from Linsley et al., 2008); SOI = Southern Oscillation Index; SPCZI = South Pacific Convergent Zone Instrumental (from Higgins et al., 2020).

Trace element analyses

Trace element analyses were carried out at the School of Earth Sciences, the University of Melbourne, Australia, using a 193 nm ArF excimer laser-ablation system coupled to an Agilent 7700 quadrupole ICPMS (Woodhead et al., 2007). Element concentrations were measured from the pre-ablated surface using a rectangular (200 × 20 μm) laser spot shape scanned parallel to the speleothem growth axis. Measurements were carried out at a speed of 20 μm/s and a laser pulse rate of 10 Hz, generating data points every ~6 μm for Pu17 and ~4.6 μm for Pu4. Five parallel scans, each spaced 5 mm apart, were performed for Pu17, while Pu4 only allowed for two adjacent parallel line scans (see Fig. S1). The repeated tracks allowed the identification of outliers caused by photomechanic ablation artifacts, which were then removed from the data. Quantification was carried out using NIST610 (Mg) and NIST612 (all other elements) glasses as external standards. Raw spectrometry data were then reduced using Iolite software and internally normalized to ⁴³Ca. The same analytical settings were used to analyze four bedrock samples collected from both inside the cave and at the surface above.

Synchrotron radiation micro-X-ray fluorescence (SR-μXRF) microscopy was performed on polished stalagmite samples at the XFM beamline at the Australian Synchrotron (Paterson et al., 2011) equipped with a Maia 384 detector array mounted 10 mm away from the sample target. The beam spot size was 1.5 μm and the monochromatic incident energy was set at 18.5 keV. The XFM spectral data were analyzed using the GeoPIXE software suite, quantified by using single element Mn, Fe and Pt foils (Micromatter, Canada), and matrix-corrected against Ca concentration (Borsato et al., 2021).

Dripwater analysis

Dripwater samples were analyzed at the University of Waikato for major and trace cations and major anions. Dissolved inorganic carbon (DIC) alkalinity was determined in situ by titration against pH using a Hach digital titrator equipped with 0.16 N H₂SO₄ cartridges. Titration data were then converted to DIC concentrations using the Gran method (Gran, 1952).

Water samples for trace elements and major cations were collected in 15 ml polypropylene tubes that were known to have low metal blanks (Nava-Fernandez et al., 2020), acidified to 2% HNO₃ (using in-house, double Teflon-distilled acid), and refrigerated until analysis. Elemental and major cation concentrations were measured using an Agilent 8900 triple quadrupole ICP-MS calibrated using external reference materials and optimized for maximum sensitivity on a daily basis, ensuring that oxides and double-charged species were less than 2%. External calibration standards were prepared using an IV71-A multi-element standard of 0.1–500 ppb for trace elements. Single-element standards were used to prepare calibration standards for major elements Ca, Fe, Si, P, S, K, and Na. An internal standard containing Sc, Ge, Te, Ir, and Rh was used for all samples. Check standards were analyzed every 20 samples and re-calibration was performed every 100 samples. Blank samples were analyzed every 10 samples to ensure minimal carryover between analyses. Repeat analyses yielded relative standard deviations (RSDs) that were < 5%.

Data processing

Several statistical and signal-processing approaches were implemented prior to establishing the transfer function between trace element concentrations and hydroclimate above the cave. The calculated infiltration was first acquired by applying the Thornthwaite equation (Thornthwaite, 1948). Potential evapotranspiration (PET) was calculated for four islands from the Southern Cooks group (Mangaia, Aitutaki, Rarotonga, and Atiu), where records of precipitation were available (Faraji et al., 2021). The calculated effective infiltration amount (P-PET) was then obtained by subtracting PET from the rainfall amount. Due to data limitations, the average temperature for each month in Atiu for the period 2014–2018 recorded outside the cave was used to calculate the entire series of PET (1914–2019 CE) and then infiltration. The infiltration values of Mangaia, Aitutaki, and Rarotonga were averaged and then combined with Atiu's record to establish a continuous calculated infiltration series from January 1914 to March 2019.

The calculated infiltration was then compared to trace elements to select the best element(s) for establishing the transfer function. Prior to making this comparison, concentrations of trace elements were detrended using wavelet decomposition. The justification for such data treatment is that each trace element shows a different general trend related to each element's specific response to aquifer-related processes. This is further supported when comparing the time series of Mg in Pu17 with that in Pu4, each showing different trends (Supplementary Information, Fig. S2). Detrending by using wavelets ensures removing only those trends that are of the lowest frequency, which are likely related to long-term aquifer process (adjustments). Wavelet decomposition detrending begins with decomposing the trace element concentration signal at nine levels using a number 4 “symlet” wavelet (Daubechies, 1992a, b). That is followed by removing level-nine (lowest frequency) information and building a new concentration signal that excludes the main trend. Since detrending shifts the baseline of the signal to zero, the baseline needs to be shifted to the previous line via adding the initial baseline value to the detrended signal. Although the detrended concentration signal is different from the original signal in terms of absolute values, each signal's sub-decadal variability is unaffected. The wavelet decomposition was implemented in MATLAB software using the Signal Multiresolution Analysis toolbox.

The detrended time series of trace elements were then resampled at a monthly resolution, allowing for better comparison between trace elements and calculated infiltration. The Pu17 time series starts in October 1672 and ends in March 2019. The Pu4 time series cover October 1800 through October 2014. After building the monthly time series, a six-month moving average based on the wet (November–April) and dry (May–October) seasons was taken. The same averaging was made on the calculated infiltration record. Since the observed precipitation, and thus the calculated infiltration for Atiu, is available only from January 1985 through March 2019, this interval was selected and used to explore the relationship between trace elements and infiltration. Through regression analysis, we selected the best trace element(s) to construct the transfer function.

RESULTS

Chronology

The Pu17 chronology is anchored by 10 U–Th ages, two of which featured large errors. Through the combined approaches presented above, this resulted in an age model with an average 2% uncertainty (Supplementary information, Fig. S1). The chronology of the Pu4 stalagmite was also constructed following the methodology described in (Faraji et al., 2021), which includes constraining the U–Th chronology using an integration of fabrics, annual hydrological cycles of LA-ICP-MS trace elements, and synchrotron μXRF two-dimensional maps of Sr. According to the final chronologies, the stalagmite Pu17 grew for 347 years with a mean growth rate of 144 ± 5 μm/year, from 1672 CE until it was collected in 2019 CE, and the stalagmite Pu4 grew for 214 years, from 1800 to 2014 CE (collection date), with a mean growth rate of 148 ± 5 μm/year, both with no apparent growth interruptions (Supplementary Information, Fig S1).

Dripwater trace elements

The drip rate and dripwater composition from five drip points in Pouatea cave were measured during the dry (October 2018) and wet (March 2019) seasons (Table 1) along with mean seawater composition (Pilson, 2012). The drip rates are higher during the wet season and the X/Ca elemental ratios are systematically higher in the dry season. Sample PUw4 is the drip point of stalagmite Pu4. The drip rate for PUw4 was 33.0 and 19.1 seconds per drip in October 2018 and March 2019, respectively. X/Ca ratios from PUw4 were always high/higher than at most other drip points. During the dry season (October 2018), PUw4 was found to have the highest Mg/Ca (0.18 mol/mol) and Ba/Ca (2.63 μmol/mol), and second highest Sr/Ca (0.89 mmol/mol). However, during the wet season (March 2019), PUw4 showed the highest Sr/Ca (0.75 mmol/mol) and second highest Mg/Ca (0.14 mol/mol) and Ba/Ca (2.30 μmol/mol) and Na/Ca (1.08 mol/mol).

Table 1. The dripwater composition from several points during the dry (October 2018) and wet (March 2019) seasons, and mean seawater composition (Pilson, 2012). The X/Ca elemental ratios are systematically higher in the dry season.

Speleothem trace elements

LA-ICP-MS trace element analyses revealed that Pu17 and Pu4 are characterized by similar concentrations in Mg, Sr, Na, and U, while P shows considerably higher values in Pu17 (Fig. 4). Pu4 is characterized by slightly higher concentrations of yttrium (Y), where more than 40% of Y concentration is below the detection limit of LA-ICP-MS analysis. Although Pu17 was fed by a slow dripping stalactite and Pu4 by a much faster one, such similarity (except for P) in the concentration of trace elements indicates that Pu17, which is the main focus of this study, is representative of the Pouatea cave. Differences in P concentration between Pu17 and Pu4 could be attributed to variations in drip rates, dissimilarities in soil type and extent above the two stalagmites at the surface, the presence of organic matter, or microbial activity.

Figure 4. Boxplots showing the ranges of trace element concentration in the two studied stalagmites. The detection limit for each element is shown as a red dash.

We compared the trace element composition spanning eight complete years of growth in the topmost section of Pu4 stalagmite to the corresponding dripwater in order to calculate the site-specific partition coefficients. The incorporation of trace elements in speleothem calcite is a multi-faceted mechanism, controlled primarily by the partition coefficient where complexation by dissolved organic material can be excluded (Hartland and Zitoun, 2018). Partitioning (D) of element X is empirically defined as: the molar ratio [X/Ca] in the calcite divided by the molar ratio [X/Ca] in the parent water (D_x =$\;{{{[ {X/Ca} ] }_{CaCO3}} \over {{[ {X/Ca} ] }_{aq}}}$) (Morse and Bender, 1990). Table 2 lists the molar ratios of Mg, Sr, Ba, and Na to Ca in the analyzed dripwater and speleothems, along with the calculated partition coefficients.

Table 2. The molar ratio of trace elements to Ca in Pouatea dripwater and speleothems, along with the calculated partition coefficient for each element. The slope of Incongruent Dissolution and Prior Calcite Precipitation vectors was calculated based on the partition coefficients of Mg and Sr using (Kd_Sr −1)/(Kd_Mg −1) (Sinclair, 2011).

(Kd_Sr −1)/(Kd_Mg −1) = 0.97 ± 0.09

Temp = 23.6 ± 0.2 °C

Moreover, the infiltrating water can experience incongruent dissolution (ICD) and/or prior calcite precipitation (PCP), shifting the dripwater composition along ICD/PCP lines, which can be characterized by a slope given by (Kd_Sr − 1)/(Kd_Mg − 1) (Sinclair, 2011). Using published values of partition coefficients of Sr (Kd_Sr) and Mg (Kd_Mg), Sinclair (2011) calculated this slope to be ~0.88 ± 0.13, which is the same as that predicted for prior calcite precipitation (PCP) (Sinclair et al., 2012). In our case, by utilizing the Pu4 partition coefficient for Mg and Sr (see Table 2), the slope of the ICD–PCP line is 0.97 ± 0.09.

Time series of trace elements were built for Pu17 and Pu4 based on the final chronologies. Time series of trace elements in Pu17 are displayed in Figure 5, along with the 12-month moving average of the calculated record of infiltration (rainfall minus potential evapotranspiration) for the years 1914–2019 (see also Supplementary Information, Fig. S3). Comparison of time series shows that in Pu17, Mg, Na, and P were best correlated with infiltration, with Mg and Na having a positive, and P having a negative, relationship. Correlation coefficients of Mg, P, and Na stand out as the best among all six elements in Pu17, being −0.64, 0.53, and −0.48, respectively, and all at statistically significant levels (p-value < 0.05) (Table 3). In Pu4, the correlation coefficient between Mg and calculated infiltration is −0.65 (p-value < 0.05) and similar to that in Pu17. However, the other elements demonstrate a poor correlation with calculated infiltration. Since Mg, Na, and P show the highest correlation with calculated infiltration in Pu17, we then computed a principal component analysis of these three elements and examined the correlation between principal component 1 (PC1), which explains 65% of the variance, and calculated infiltration. The regression analysis of the six-month moving averaged PC1 and calculated infiltration returned a coefficient of −0.60, which is also statistically significant.

Figure 5. Comparison of time series of trace elements (Mg, Sr, Na, U, P, Y, and PC1 of Mg–Na–P PCA) in Pu17 with the 12-month moving average of the calculated infiltration (rainfall minus potential evapotranspiration) from 1914 to 2019. Mg, Na, and P have the best correlation with infiltration, with Mg and Na showing a negative relationship, and P showing a positive relationship (see Table 3). It seems that low-infiltration years are less well captured by trace elements. Please note that the y-axes are in reversed order for Mg, Na, U, and PC1. The red triangles next to the y-axes show the decreasing or increasing order of the y-axes.

Table 3. The correlation coefficients and p-values for the least-squared regressions between each trace element and calculated infiltration for all the studied elements in Pu17 and Pu4. PCA Mg-Na-P PC1 is a principal component analysis of these three elements to examine the correlation between principal component 1 (PC1), which explains 65% of the variance, and calculated infiltration. The correlation coefficients better than ± 0.5 are highlighted in bold.

DISCUSSION

Sources of trace elements in Pouatea dripwater

Dripwater has been studied extensively in inland caves, providing a basis for interpreting the speleothem trace element record. This has been achieved by constraining the origin of trace elements in cave dripwater through the chemical characterization of endmembers (bedrock, soil, dust, etc.) (e.g., Goede et al., 1998; Musgrove and Banner, 2004; Tremaine and Froelich, 2013). However, Pouatea is a coastal cave. Unlike inland caves far from the coast, sources of ions, especially for Na and Mg, at coastal caves may also include sea-spray aerosols delivered by meteoric deposition to the epikarst. In their study of a coastal cave on Niue island, tropical South Pacific, Tremaine et al. (2016) discussed and quantified the marine aerosol contribution for different elements, showing that sea-salt input accounted for a large portion of dripwater Na and Mg (89% and 85%, respectively) and a smaller portion of the dripwater Ca and Sr (19% and 17%, respectively). We, therefore, anticipated a significant contribution from marine aerosols to Pouatea dripwater, given its proximity (0.5 km) to the coast. The possible provenances of different trace elements in Pouatea cave dripwater were briefly evaluated through a Principal Component Analysis and discussed in Faraji et al. (2021), which suggested that most Sr originates from the host rock and that Mg and Na are related to marine aerosols, whereas Y appears to be a soil-derived element similar to P (Hartland et al., 2012). The origin of P concentrations in tropical South Pacific islands appears to be related to the weathering of corals (Whitehead et al., 1993). It exhibits a complex behavior, possibly resulting from the influence of microbial activity at the surface of a growing speleothem (Frisia et al., 2012).

In order to better understand the potential sources of trace elements in Pouatea dripwater and investigate the trace element incorporation in Pu4 and Pu17 stalagmites, we developed a conceptual model following the approach outlined in Sinclair (2011). At first, we measured the dripwater composition from several points during the dry and wet seasons and compared these to the seawater and host rock Mg/Ca versus Sr/Ca of the several host-rock samples that we collected at the surface (Fig. 6). It is essential to note that the Pouatea host rock is very porous and partially diagenetically altered, and therefore its composition can vary substantially from site to site. Consequently, the host rock sampled in the cave near the stalagmite is unlikely represent the actual composition of the epikarst where most of the host rock dissolution occurs. For this reason, it is preferable to analyze several rock samples above the cave (HR1–HR4) in order to document the range in the trace element composition of the host rock in the infiltration area. HR1 was sampled from inside Pouatea cave, HR2, HR3, and HR4 were sampled from above the Pouatea cave.

Figure 6. The natural logarithm of the molar ratio of Mg/Ca versus Sr/Ca for the analyzed host rocks (HR), dripwater, mean seawater, and the top 2 mm of stalagmites Pu17 and Pu4. Dripwater corresponding to stalagmite Pu4 (PUw4) is plotted as light blue circles. The two blue curves (seawater mixing curves) show the composition of dripwater evolved from the congruent dissolution of HR2 and HR3 host rock with an increasing amount of seawater contribution. Dots on the seawater-mixing curves show the incremental (0.1%) contribution from sea spray. The slopes of the seawater-mixing curves for HR2 and HR3 in the seawater range from 0.2% to 0.7% and vary between 0.22 and 0.44. The red line shows the dripwater evolution following PCP starting from a HR2 composition and 0.2% seawater contribution. Dots on the PCP line show the incremental (10%) occurrence of PCP (to 70%), which gradually enhances the ln (Mg/Ca) and ln (Sr/Ca) of dripwater along a 0.97 slope. The slopes (s.) for each pair (dry–wet season) of dripwater analyses are reported in the legend. The slopes of PUw2 and PUw3 are aligned with the seawater vector, PUw4 and PUw6 have slopes intermediate between the seawater and PCP vectors. PUw5 has a slope exceeding the slope of the PCP vector, suggesting contribution from different HR during the wet and the dry seasons. Therefore, dripwater in Pouatea cave undergoes a complex evolution governed by different amounts of sea-spray contribution after dissolution of host rock of different composition, and, finally, a variable amount of PCP. The entire dripwater dataset can be reproduced with HR composition ranging from HR2 and HR3, sea-spray contribution 0.2–0.7%, and PCP 0–30%.

Figure 6 shows the cross plot of the natural logarithm of elemental molar ratios of Mg/Ca versus Sr/Ca for seawater, host rocks (HR1–HR4), dripwater, and the topmost 2 mm of the Pu17 and Pu4 stalagmites. Compositions of the four different host rocks vary substantially, particularly with respect to Sr/Ca ratios (Fig. 6). This means that dissolution in different parts of the epikarst is expected to produce dripwater with a wide range of initial Sr/Ca, and a limited range in Mg/Ca. The comparison between host rock and Pouatea cave dripwater indicates that dripwater is generally characterized by similar Sr/Ca ratio, but a much higher Mg/Ca ratio, suggesting a substantial seawater contribution. Hence, we modeled the dripwater composition starting from different host rock compositions and mixing an increasing amount of seawater. The simulation suggests that Pouatea cave dripwater can be reproduced by the congruent dissolution of host rock, with a composition between HR2 and HR3, and a limited seawater contribution (0.2–0.7%). In particular, dripwater PUw4 (feeding the Pu4 stalagmite) can be modeled with the congruent dissolution of host rock HR3 and 0.53% seawater contribution, whereas dripwater from other cave sites possibly reflects a host rock composition intermediate between HR2 and HR3. Interestingly, the host rock HR1 sampled inside the cave near Pu4 and Pu17 stalagmites has the lowest Sr/Ca and Mg/Ca ratios and is not representative of the dripwater composition.

Moreover, infiltrating water can be modified further by incongruent carbonate dissolution (ICD) and/or prior calcite precipitation (PCP), thus shifting the dripwater composition along ICD/PCP lines characterized by a slope given by (Kd_Sr − 1)/(Kd_Mg − 1) (Sinclair, 2011). By plotting the ICD/PCP lines in ln (Mg/Ca) versus ln (Sr/Ca) space, we observe that Pouatea dripwater can be simulated by using a host rock composition intermediate between HR2 and HR3 and variable amounts (0–30%) of PCP. However, the almost constant Na/Ca composition in dripwater suggests a narrow range in seawater contribution of 0.41–0.71% (Table 4) and, therefore, a limited contribution of ICD–PCP.

Table 4. Host rock, mean seawater and PUw4 dripwater composition. The dripwater modeling show the theoretical dripwater composition derived from the dissolution of a host rock with intermediate composition between HR2 and HR3. The host rock and seawater contributions to PUw4 dripwater were calculated based on Ca and Na content.

¹(Pilson, 2012)

²Dripwater emerging from the dissolution of a host rock with intermediate composition between HR2 and HR3.

³Seawater contribution of 0.533% accounts for 2.563 mmol of Na in PUw4 dripwater, while host rock dissolution contributed only 0.0007 mmol of it.

In order to explain the dripwater trace element composition we, therefore, needed to take into account the contribution of marine aerosols that shift the dripwater composition towards mean seawater. For this purpose, we used Na⁺ as a conservative seawater tracer (see Tremaine et al., 2016). By assuming that almost the entire amount of Na⁺ in dripwater originates from marine aerosol, we calculated the percentage of sea-spray contribution to the trace element composition of dripwater (see Table 4). Our model suggests that 90% of Mg, 30% of Sr, and less than 5% of Ba in the Pouatea dripwater originate from sea spray, while 99% of Ca, 10% of Mg, 75% of Sr, and 41% of Ba are sourced from the host rock. The sum of seawater and bedrock contribution accounts for nearly 100% of the Ca, Ma, Na, and Sr, but only accounts for 45% of Ba. It has been well documented that the concentration of Ba in sodic Oceanic Island Basalts (like Atiu) is extremely high compared to all other oceanic basalts (Wood, 1978; Niu and O'Hara, 2003). Therefore, leachates from the weathered volcanic core of the Atiu island likely are the other source that contributes to the supply of Ba into cave dripwater.

The two blue lines on Figure 6 show the evolution of dripwater derived from the congruent dissolution of host rock samples and an increasing sea-spray contribution. The slopes of the seawater-mixing lines vary from 0.20 to 0.45, depending by the initial HR composition and by the amount of sweater contribution, while the slope of the PCP lines, which is based on the partition coefficient of Mg and Sr calculated from PUw4 − PU4 (see Table 2), is 0.97. When water analyses from the same drip site are aligned along a seawater-mixing line (as the case of PUw2 and PUw3 dripwater) they suggest the congruent dissolution of a host rock with homogeneous composition, a variable contribution of marine aerosols, and no significant modification by PCP. On the other hand, when water analyses from the same drip site diverge from the seawater curves with a slope approaching 0.97, they suggest an increasing contribution of PCP (PUw4, PUw6), while dripwater slopes exceeding 0.97 (PUw5) can be attributed to the dissolution of host rocks with different compositions during different seasons (i.e., under different saturation states of the aquifer).

Incorporation of trace elements in the stalagmites

Most studies focus on elements that form divalent cations in solution and substitute for Ca in the carbonate crystal lattice, particularly Mg, Sr, and Ba (together with Mn and Fe in reducing waters). The partition coefficients of these elements may vary to a greater or lesser extent with temperature, precipitation rate, crystal morphology, or other aspects of solution composition (Fairchild and Treble, 2009). Other elements may be incorporated in calcite by different mechanisms. For example, Borsato et al. (2007) documented that a variety of trace elements (Y³⁺, Zn²⁺, Cu²⁺, Pb²⁺, P³⁺, Br⁻) enrichment coincides with visible laminae containing fluorescent organic matter. Also, soil bacteria may enhance the mobility of certain trace elements (Fe²⁺, Cu²⁺, Ni²⁺, Co²⁺, Mn²⁺, Y³⁺) through the karst system by releasing chelating compounds (Hartland et al., 2012), while at the same time inhibiting the partitioning of the trace elements into stalagmites (Hartland et al., 2014; Hartland and Zitoun, 2018). However, because of the interactions of trace elements, particulate matter, and organic colloids, one simple model cannot be assumed to explain trace element partitioning independently onto their preferred crystal sites (Hartland et al., 2012; Pearson et al., 2020). Therefore, connecting trace element distributions with organic-stimulated fluorescence requires care (see Sliwinski and Stoll, 2021).

A number of laboratory and cave studies have demonstrated a temperature dependence on D_Mg (Gascoyne, 1983; Huang and Fairchild, 2001; Day and Henderson, 2013; Wassenburg et al., 2020). Sliwinski and Stoll (2021) reported a lateral variability in Mg distribution associated with fluorescent features. They hypothesized that although Mg partitioning in calcite is controlled by site preference and inhibition by colloids, low concentrations of dissolved organic matter (DOM) and high supersaturation (as in Pouatea dripwater) lead to uniform organic occlusion and uniform trace element distribution. This is further supported by recent experimental findings showing that DOM adsorption heavily influences the crystal structures with finer, smooth-faced, rhombohedral crystals that form in growth solutions with low aqueous DOM, and prismatic, ‘impure’ crystals produced at high aqueous DOM (Pearson et al., 2020). Others have reported growth rate dependencies on D_Na, D_Sr, and D_Ba in both laboratory calcite and in speleothems (Tesoriero and Pankow, 1996; Huang and Fairchild, 2001; Treble et al., 2003, 2005b; Gabitov and Watson, 2006; Day and Henderson, 2013; Belli et al., 2017). Wassenburg et al. (2020) observed a positive correlation between D_Sr and the calcite Mg/Ca ratio, confirming the earlier findings from Mucci and Morse (1983). Additionally, recent studies (e.g., Drysdale et al., 2019) suggest that the ratio of Na/Sr (mmol/mmol) in the parent water controls the partition coefficient of Sr. Drysdale et al. (2019) documented an exponential decrease in D_Sr with increasing Na/Sr described by the function: D_Sr = 0.1173e^{−0.0004(Na/Sr)}. Given the high concentration of Na in Pouatea dripwater, changes in the Na content will directly affect the dripwater Na/Sr ratio, in turn, modifying D_Sr, and hence the slope of the precipitation line (dissolution-PCP-precipitation line; Sinclair, 2011; Sinclair et al., 2012). For this reason, the modeling of the Mg/Ca versus Sr/Ca system has to take into account variability of the Na/Sr ratio in the dripwater.

To examine the trace element composition of Pu17 and Pu4 stalagmites, three points were selected along the seawater curves (Fig. 7) that represent solutions with 0.2%, 0.4%, and 0.7% sea-spray contribution, respectively. Then, we calculated the theoretical composition of precipitated calcite based on the in-situ partition coefficients (Table 4). The amount of sea-spray contribution directly controls the Na/Sr dripwater ratio and, subsequently, the partition coefficient of Sr (Drysdale et al., 2019). This, in turn, modifies the slopes of the ln (Mg/Ca) versus ln (Sr/Ca) precipitation lines, which vary from 0.70 to 0.85 (corresponding to a D_Sr of 0.050 and 0.084, respectively), depending on the amount of sea-spray contribution. The endpoints of the three precipitation lines account for the entire Mg/Ca range in the stalagmites, while different host rock sources and/or additional ICD/PCP are required in order to account for the entire range in Sr/Ca (Fig. 7).

Figure 7. The natural logarithm of molar ratio of Mg/Ca vs. Sr/Ca for HR2 and HR3 host rock, PUw4 dripwater, and the entire time series of Pu17 and Pu4 stalagmites. The two blue curves (seawater mixing curves) show the composition of dripwater evolved from the congruent dissolution of HR2–3 and HR3 host rocks, with an increasing amount of seawater contribution. Dots on the seawater-mixing curves show the incremental (0.1%) contribution from sea spray. In order to model the stalagmite composition, we selected three points corresponding to the dissolution of HR2–3 host rock with 0.2%, 0.4%, and 0.7% sea spray contribution, respectively, and constructed the precipitation lines on the basis of the calculated D_Mg and D_Sr (Table 2). The calculated D_Sr is influenced by the dripwater Na/Sr ratio and varies from 0.050 (seawater contribution 0.7%) to 0.084 (seawater contribution 0.2%). This, in turn, is reflected by changes in the slope of the precipitation lines that vary from 0.70 to 0.85 as a function of increasing sea-spray contribution. As shown, the three selected seawater contributions can reconstruct the composition of the entire datasets in the two stalagmites.

Moreover, the correlation coefficients for the entire Mg/Ca versus Sr/Ca datasets of the two stalagmites are not significant (r² = 0.01 and 0.18 for Pu17 and Pu4, respectively), and their slopes are more aligned with the seawater slopes (0.2–0.45) and not with the ICD/PCP slope (~0.8). This indicates that seawater contribution is the main controlling mechanism for Mg incorporation in the stalagmites, while Sr incorporation is also influenced by ICD/PCP. Then, given the sources of the other analyzed elements (listed in Table 4), we confidently can affirm that the incorporation of Mg and Na in the stalagmites is primarily modulated by differential sea-spray contribution, while Sr and Ba are mainly controlled by ICD/PCP.

Finally, the sea-spray contribution is diluted by the rainfall amount and, therefore, a higher Na and Mg content in dripwater/stalagmites indicates dryer conditions. At the same time, drier conditions enhance the ICD/PCP contribution, and hence the Sr and Ba content. Therefore, as pointed out in Faraji et al. (2021), rainfall and infiltration amounts are the ultimate mechanisms controlling the concentration of trace elements in Pouatea speleothems at a seasonal scale, regardless of the trace element provenance. Higher infiltration (i.e., higher rainfall at a relatively constant temperature) corresponds to less rock-water interaction (RWI), less PCP, and dilution of marine aerosol-derived elements in dripwater. Thus, when the drip rate is relatively high, formation of calcite layers with a more porous fabric and relatively low Mg, Na, Sr, Ba concentrations is expected. By contrast, in wet periods, soil-derived elements, such as Y (and, in part, P), would increase in the absence of other process controls, such as organic complexation (Hartland and Zitoun, 2018; Sliwinski and Stoll, 2021). It is, therefore, expected for infiltration to hold a negative correlation with Mg, Na, Sr, Ba and a positive correlation with Y and P. However, among these elements, Sr is the most complicated and is controlled by four different factors: growth rate, PCP, dripwater Na/Sr ratio, and calcite Mg/Ca ratio. The interplay of these factors modulates the concentration of Sr in speleothem calcite, which, depending on the significance of each factor, can show either a positive or negative correlation with infiltration. For example, if the growth rate is the governing mechanism, then Sr would show higher values during wet seasons (positive correlation with infiltration), whereas if PCP and calcite Mg/Ca ratio are the controlling factors, the correlation with infiltration would be negative. In the case of Pu17 (see Fig. 5 and Table 3), Sr is positively correlated with the infiltration, which suggests that growth rate can outcompete the role of PCP.

A trace element record of infiltration for the last ca 350 years

As previously discussed, Mg held the highest correlation with infiltration, and therefore, was selected as a proxy to reconstruct infiltration in Atiu. Mg stands out as the best element for the reconstruction of infiltration because it originates mostly from a single source (marine aerosols). In contrast, Sr has a double provenance (host-rock dissolution and, to a minuscule extent, marine aerosols). Plus, Sr incorporation in stalagmite calcite is controlled by contrasting mechanisms such as growth rate, PCP, Na/Sr ratio, and calcite Mg/Ca ratio. Additionally, the composition of different host rocks inside (HR1) and above (HR2, HR3, and HR4) Pouatea cave revealed that solutions emerging from these host rocks are expected to have substantially different Sr/Ca, but fairly similar Mg/Ca ratios, which is important when considering the dripwater in terms of dry and wet seasons. During the dry season, a narrow area of host rock will be affected by infiltration, whereas during the wet season, a more extensive degree of water–rock interaction would generally occur. Given the difference between Mg/Ca and Sr/Ca ratios of the different host rocks, such a contrast between dry and wet can substantially alter the dripwater Sr/Ca ratio, while the Mg/Ca ratio is dominated by the amount of sea spray. Another reason that Mg proved to be a better proxy than Sr, is that the dripwater Na/Sr (Drysdale et al., 2019) and calcite Mg/Ca ratio (Wassenburg et al., 2020) control the partition coefficient of Sr (D_sr), whereas D_Mg is indifferent to the dripwater Na/Sr ratio.

Magnesium was also superior to Na as a proxy for paleohydrology, even though, like Mg, Na mainly derives from a single source. This disparity can be attributed to the incorporation of Na⁺ in calcite, which mostly depends on the availability of intercrystalline sites (Busenberg and Plummer, 1985), whereas Mg²⁺ substitutes directly for Ca²⁺ in the calcite lattice (Mucci and Morse, 1983) and is mostly unaffected by crystal fabric. Therefore, Mg is anticipated to be the best hydrologically sensitive element. Figure 8 illustrates the transfer function between six-month averaged Mg and calculated infiltration for both Pu17 and Pu4. Although the two stalagmites have different characteristics in terms of drip rate, growth rate, and petrography, the slopes of the correlation between Mg content and calculated infiltration are similar. This finding is fundamental to this study because it shows that Pu17 and Pu4 responded similarly to the climate drivers despite several differences in their characteristics.

Figure 8. The transfer function between six-month averaged Mg and calculated infiltration for both Pu17 and Pu4. Although the two stalagmites have different characteristics in terms of drip rate, growth rate, and petrography, they show a similar slope in the transfer function. This consistency is fundamental for this study because it shows that Pu17 and Pu4 responded similarly to the climate drivers above the cave despite some differences in their characteristics. The Pu17 function was used to reconstruct infiltration from the time series of Mg in Pu17. The color scale shows the dry (red) and wet seasons (blue).

Having established the transfer function between the concentration of Mg and infiltration for the two stalagmites, we then applied the transfer functions on the time series of Mg in Pu17 and Pu4 to reconstruct the variability of infiltration. The top panel of Figure 9 compares the calculated infiltration (12 months average) with reconstructed infiltration for Pu17 and Pu4 (see also Supplementary information, Fig. S4). Despite differences between reconstructed and calculated records in absolute values, it is clear that the variability of infiltration reconstructed using Pu4 and Pu17 Mg is similar to the variability observed in the infiltration record calculated from instrumental data. For most parts of the record, the transfer function is typically able to reconstruct relatively dry periods accurately. The discrepancy between reconstructed and calculated values for the older part of the record (1914–1930 CE) can be attributed to the uncertainties associated with the transfer function in that part. However, it should be noted that the instrumental rainfall records (and calculated infiltration) up to 1930 are derived solely from Mangaia, and as such, the calculated infiltration for Atiu is more uncertain during this period.

Figure 9. Reconstructed infiltration for Pu17 and the instrumental records of SOI, Niño 3.4 SST, SPCZI, IPO, and TCs (tropical cyclones) count. The figure shows that dry events recorded in the instrumental climate data are also traceable in reconstructed infiltration, even though the strength of the events may vary between different records. Intervals highlighted in yellow on the bottom panel indicate dry periods where our reconstruction overestimates the strength of the events, which likely is caused by TCs that supply extra Mg to the cave system. The red-outlined triangles are tropical cyclones that likely caused extreme rainfall events, followed by dilution of trace elements in the aquifer. Storms and TCs affecting Southern Cooks [1] data from de Scally, 2008; d'Aubert and Nunn, 2012; Diamond et al., 2012.

In Figure 9, the reconstructed infiltration derived from Pu17 Mg series is also plotted against the well-known drivers of rainfall in the South Pacific (SOI, Niño 3.4 SST, SPCZI, and IPO) for the period 1870–2019. During positive-phase IPO/El Niño, the SPCZ displaces northeastward (zonal) and the Southern Cook Islands (SCI) are relatively dry. Other phases (positive-phase IPO/La Niña, negative-phase IPO/El Niño, negative-phase IPO/La Niña) mark wet conditions in SCI corresponding to a southwestward (diagonal) SPCZ. Figure 9 also includes the infiltration reconstructed using P and U transfer functions and years when the Southern Cook Islands were affected by storms and/or tropical cyclones (TCs) (data from de Scally, 2008; d'Aubert and Nunn, 2012; Diamond et al., 2012). The reconstructed infiltration detects most dry periods recorded based on SOI and Niño 3.4 SST. For example, the strong El Niño events (dry events) of 1886–1888, 1912–1914, 1940–1942, 1965–1966, 1982–1983, and 1997–1998 appear on the reconstructed record as decreased infiltration (Fig. 9).

The relationship between SPCZI and reconstructed infiltration might seem less clear. Unfortunately, very little can be inferred from the magnitude of the SPCZI because the SPCZI varies monotonically but non-linearly with displacement of the SPCZ, especially in the region southeast of Samoa where the Southern Cook Islands are located (Weir et al., 2021). There are examples where reconstructed infiltration overestimated the severity of dry periods (e.g., the dry event of 1967–1970). This could be attributed to furious winds and storms reported for the Southern Cook Islands in 1967 and 1970 (d'Aubert and Nunn, 2012) that can supply extra Mg to the cave system, which eventually can appear as exaggerated dry periods. Interestingly, the infiltration reconstructed from P and U, which are unaffected by sea spray, show higher values of infiltration than the Mg record for the period 1967–1970. Intervals highlighted in yellow on the bottom panel in Figure 9 point to dry periods where our reconstruction overestimates the strength of the events, likely caused by TCs that supply extra Mg to the cave system (e.g., 1913, 1941–1944, 1953–1954, and 1962–1968).

Tropical cyclones can also cause extreme rainfall events, which inevitably are followed by dilution of trace elements in the aquifer. These events appear as wet periods on the Mg, P, and U reconstructed records of infiltration (e.g., 1905, 1956, and 1959 storm/cyclone years in Fig. 9, displayed as red-outlined triangles). The 10-year moving average of reconstructed infiltration seems to follow the IPO phases for the period 1870–2016, which is consistent with the decadal effect of IPO in modulating rainfall in the Southern Cook Islands (Fig. 9).

To further evaluate the accuracy of the reconstructed infiltration, we draw on available reconstructions of the ENSO anomaly (Li et al., 2011), SPCZI (SPCZI_r) established by Higgins et al. (2020), and IPO (Linsley et al., 2008; Vance et al., 2022). These records are shown and compared with reconstructed infiltration from Pu17 in Figure 10. Interestingly, our reconstructed infiltration matches reasonably well with these records, even though discrepancies exist because these are reconstructions and, undoubtedly, bear uncertainties. However, clear agreements between reconstructed infiltration and ENSO (Fig. 10a) and SPCZI_r (Fig. 10b) support the accuracy of our reconstruction. This is of great importance because ENSO drives much of the interannual variation in rainfall (and thus infiltration) in the Southern Cook Islands (Weir et al., 2021) through displacement of the SPCZ, with the SPCZ moving northeastward during El Niño events and southwestward during La Niña events (Salinger et al., 2014). Therefore, the variability of infiltration reconstructed in this study can be interpreted as rainfall variation caused by the movements of the SPCZ following ENSO events. The 10-year moving average of reconstructed infiltration (Fig. 10c) shows variabilities very similar to the IPO reconstructed by Linsley et al. (2008) and mostly different from the record by Vance et al. (2022). In addition to the ENSO-triggered displacement of SPCZ, there is also a long-term movement associated with the IPO, with the SPCZ moving to the northeast in the positive phase and to the southwest in the negative phase of IPO (Folland et al., 2002; Salinger et al., 2014). The similarity between reconstructed infiltration at decadal scales and IPO might indicate the decadal influence of IPO on SPCZ that leads to long-term rainfall variabilities.

Figure 10. Infiltration derived from the time series of Mg in Pu17 compared with ENSO anomaly (Li et al., 2011), reconstructed SPCZI (SPCZI_r) (Higgins et al., 2020), and IPO (Linsley et al., 2008; Vance et al., 2022) as the main factors driving climate in the South Pacific. The figure shows that the reconstructed infiltration matches very well with these records, even though discrepancies exist because these are reconstructions and, undoubtedly, bear uncertainties.

Therefore, the reconstructed infiltration for Atiuan stalagmites provides vital information about hydroclimate variability in the past 350 years and paves the way to reconstruction of hydroclimate variability beyond the past 350 years, for longer periods in the past, possibly beyond existing ENSO and SPCZ reconstructions. It also open the door to an enhanced understanding of the response of the SPCZ to a globally warming climate. This study emphasizes the significance and substantial potential of trace elements in paleoclimate interpretations of speleothems.

CONCLUSION

This study presents a high-resolution reconstruction of effective infiltration in the Southern Cook Islands. We studied two active stalagmites, which both grew for the past 200 years. Stalagmite Pu17 spanning from 1672 to 2019 CE, and stalagmite Pu4, spanning from 1800 to 2014 CE, retrieved from Pouatea cave on Atiu Island. These two modern stalagmites provide an excellent opportunity to examine the potential of Atiuan speleothems as sources of information about paleoclimate and hydroclimate variability in the remote and climatically vulnerable South Pacific islands.

We used LA-ICP-MS trace elements as proxies to reconstruct the infiltration in Atiu. In order to better understand the provenance of trace elements in Pouatea cave dripwater and stalagmites, we developed a model that takes into account the contribution of trace elements from different sources and how trace elements are subsequently incorporated in speleothems. Our modeling indicated that almost all the Na and nearly 90% of Mg in Pouatea dripwater is sourced from sea spray. Sr and Ba originated mainly from bedrock, with additional Ba supplied by marine aerosols and leaching of Ba-rich weathered sodic oceanic basalt. Based on this process understanding, we determined that infiltration is the driving mechanism of most trace elements in Pouatea at a seasonal scale, regardless of provenance. Higher infiltration (i.e., higher rainfall at a relatively constant temperature) corresponds to less rock–water interaction (RWI), less PCP, and reduced marine aerosol contributions. Based on this, we established transfer functions between trace elements and calculated effective infiltration, from which Mg had the best correlation leading to its use in reconstructing infiltration.

Reconstructions based on Mg were compared with records of climate drivers in the South Pacific. Reconstructed infiltration captures dry periods in instrumental records. The good correlation between our reconstruction and the SPCZI and ENSO confirms the accuracy of our reconstruction. It also emphasizes the substantial role of these two climate drivers in controlling precipitation in the South Pacific.

Overall, our results show an unprecedented high-resolution reconstruction of infiltration based on speleothem trace elements in the South Pacific. Evaluating the results of this study against instrumental data demonstrates the potential of Atiuan speleothems to reconstruct hydroclimate and deliver well-dated and accurate records of past climate to shed light on the hydroclimate variability in the South Pacific.