Diatom–environment relationships and limnological variability: an updated quantitative tool for palaeoclimatology on sub-Antarctic Macquarie Island

Sub-Antarctic Macquarie Island is ideally located for reconstructing past variations in Southern Hemisphere westerly wind strength. Diatoms are a valuable palaeolimnological tool on sub-Antarctic islands, providing a means to reconstruct past climate and environmental changes. Diatom communities are sensitive to changes in lake electrical conductivity (EC) linked to westerly wind–driven sea-spray inputs on Macquarie Island, and diatom–conductivity models have previously been used to infer past westerly wind variability. Here we present new diatom data from 52 lakes to assess diatom–environment relationships and develop an updated diatom–conductivity model for Macquarie Island. Seasonal and multi-year water chemistry and isotope data were analysed to assess temporal variability in hydrochemical processes and the influence of evaporation, ensuring the resulting diatom-conductivity model reflects external climatic drivers rather than local dynamics. Statistically robust transfer functions were developed for EC (bootstrapped r² = 0.80, RMSEP = 0.40), while pH had weaker predictive performance. For EC, weighted averaging and maximum-likelihood approaches performed comparably, although the former showed reduced predictive power at high EC where low species turnover and nutrient collinearity affected accuracy. This quantitative-diatom model combined with understanding of hydrogeochemical processes provides an improved basis for reconstructing past Southern Hemisphere westerly wind variability, which can be applied in future palaeoclimate studies on Macquarie Island.

1 Introduction

The Southern Ocean region exerts a strong influence on Southern Hemisphere and global climates (Jones et al., 2016; Fogt and Marshall 2020). Sub-Antarctic islands are among the few landmasses located in the Southern Ocean, making them important sites for understanding the past and future role of the Southern Ocean on climate variability. The Southern Hemisphere westerly winds (SHW) are a major driver of Southern Hemisphere mid- to high-latitude climates, modulating ocean circulation, mid-latitude temperature and precipitation regimes, and the efficiency of the Southern Ocean carbon sink (Gillett et al., 2006; Le Quéré et al., 2009; Fletcher et al., 2021; Menviel et al., 2023; Thomas et al., 2025). Instrumental data show that in recent decades the SHW have intensified and shifted poleward in response to warming (Marshall, 2003; Fogt and Marshall, 2020). These changes have been linked to an increase in net outgassing of carbon dioxide (CO₂) from deep-storage reservoirs in the Southern Ocean, with significant implications for future atmospheric CO₂ levels and global temperatures (Goyal et al., 2021; Nicholson et al., 2022; Mongwe et al., 2024; Olivier and Haumann, 2025). Understanding long-term SHW variability is key to assessing the impacts of SHW dynamics under future climate warming scenarios.

Diatoms are highly sensitive to environmental changes and are widely used as palaeolimnological proxies to infer climate and environmental changes (Roberts et al., 2000; Verleyen et al., 2004; Sterken et al., 2008; Recasens et al., 2015; Liao et al., 2020; Peng et al., 2022; Deng et al., 2025). Diatoms are well established as indicators of salinity and ionic composition, forming the basis of numerous diatom–salinity or -conductivity transfer functions across a range of environments (Gasse et al., 1997; Verleyen et al., 2003; Volik et al., 2017; Maslennikova, 2020; Farqan et al., 2025). These approaches have been successfully applied in diverse settings demonstrating the reliability of diatom assemblages for reconstructing past hydrochemical and environmental change.

Previous work on sub-Antarctic Islands has demonstrated that aquatic diatom communities are significantly influenced by changes in salinity (inferred from electrical conductivity, EC), allowing quantitative diatom-conductivity models to be developed (Gremmen et al., 2007; Saunders et al., 2009, 2015, 2018; Perren et al., 2020, 2025; Van Nieuwenhuyze, 2020). On sub-Antarctic islands, lake water salinity changes are largely controlled by wind-driven sea spray aerosol (SSA) inputs, via both wet and dry deposition, with increased inputs occurring when winds are stronger and vice versa (Evans, 1970; Buckney and Tyler, 1974; Saunders et al., 2009, 2015; Humphries et al., 2021). Based on this, diatom-conductivity transfer functions have been used to infer past Holocene SHW intensity on Macquarie Island (Saunders et al., 2018), Marion Island (Perren et al., 2020) and in southern South America (Perren et al., 2025). These relationships are understood to reflect longer-term, integrated hydrogeochemical and ecological responses to persistent wind-driven sea-spray inputs, rather than event-scale meteorological forcing.

Earlier studies on Macquarie Island have analysed diatom-environment relationships (McBride, 2009; Saunders et al., 2009) and their application as palaeoenvironmental and climate proxies (Keenan, 1995; Saunders et al., 2013, 2018; Deng et al., 2025). However, from the late 1900s to early 2000s overgrazing from increasing invasive rabbit populations (up to 150 000 individuals estimated from 2005–2006) resulted in widespread ecosystem degradation, including erosion, vegetation loss, and altered organic inputs into lakes (Scott and Kirkpatrick, 2008; Terauds, 2009). This affected aquatic ecosystems and diatom diversity (Marchant et al., 2011; Saunders et al., 2013). The Macquarie Island Pest Eradication Programme successfully eradicated all invasive vertebrates (principally rabbits) from the island in 2011, triggering substantial ecosystem recovery (Springer, 2018; Fitzgerald et al., 2021). Although direct limnological data to assess ecosystem recovery of individual lakes is not available, widespread vegetation recovery following early efforts of the eradication programme in 2010–2011 provides strong evidence that ecosystem processes across the island are no longer characterised by extreme disturbance (Shaw et al., 2011; Springer, 2018; Fitzgerald et al., 2021). This is expected to have decreased catchment erosion and sediment and nutrient delivery into lakes relative to the peak disturbance period, resulting in post-eradication (recovering) ecosystems.

Reassessing diatom–environment relationships under current post-eradication conditions is necessary, because earlier studies were conducted during a period of vertebrate-induced disturbance rather than under near-natural conditions (Saunders et al., 2013). Developing new diatom models based on post-eradication conditions may better represent pre-invasion baseline communities, improving the accuracy and ecological relevance of palaeolimnological reconstructions. Furthermore, incorporating revised taxonomy and newly identified species will enhance the model's ecological resolution and predictive performance.

Understanding the processes that drive lake water chemistry, such as precipitation, evaporation, groundwater inputs, and nutrient cycling, and how they vary across temporal and spatial scales is essential when interpreting diatom–environmental relationships. Meredith et al. (2022) defined lake hydrogeochemical processes across Macquarie Island showing that dominant processes vary locally, and lakes can be classified as predominantly influenced by SSAs, catchment processes (i.e., with greater water-rock interaction), or precipitation (i.e., more dilute lake waters). Sea-spray-influenced lakes occur near the west coast and on the western edge on the Macquarie Island plateau, where exposure to the SHW is greatest. In contrast, catchment-influenced lakes with higher terrestrial ion concentrations are found at lower elevations, and rainfall-influenced lakes with low ion concentrations occur at higher elevations. This hydrogeochemical framework supports the hypothesis that for lakes near the west coast, including those on the western edge of the plateau, EC-related diatom variation on Macquarie Island primarily reflects SHW-driven sea-spray inputs rather than local hydrological or geochemical controls.

While present day water chemistry provides valuable insight into spatial variability, it is necessary to quantify temporal variability in hydrogeochemical processes, particularly evaporation, to assess how seasonal, interannual, and longer-term changes modify ion concentrations in lakes, including those derived from SSA. Establishing seasonal and multi-year lake water hydrogeochemical datasets will enhance confidence in proxy interpretations and form a foundation for long-term monitoring of Macquarie Island lakes. Such research is rarely applied to develop ecological transfer functions particularly in such remote, isolated settings, and has not yet been undertaken on other sub-Antarctic islands. These factors highlight the importance of this work for understanding how sub-Antarctic Island ecosystems will respond to future climate and environmental changes, particularly given the rapid ecological shifts in response to climate that are already documented across the region (Le Roux and McGeoch, 2008; Lee and Chown, 2016; Nel et al., 2023).

Here, we present new data from lakes on Macquarie Island quantifying post-pest eradication relationships between surface-sediment diatom communities and environmental conditions. Using comprehensive water chemistry datasets from 2018 and 2022–2023, we examine seasonal and interannual variability to develop updated diatom–environment transfer functions. This integrated approach strengthens the application of diatom-based proxies and provides a first step towards long-term monitoring of sub-Antarctic lake systems by contributing to baseline data and establishing an analytical framework to track ecological and biogeochemical change. These transfer functions will be applied in future studies to reconstruct past Holocene climate variability on Macquarie Island.

2 Methods and Materials

2.1 Study area: Macquarie Island

Macquarie Island (54°50^′ S, 158°85^′ E) is a small sub-Antarctic island (128 km²) located in the Southern Ocean just north of the polar front, 1200 km south-west of New Zealand and 1300 km from the Antarctic continent (Fig. 1). It is one of the few landmasses within the Polar Frontal Zone and modern core SHW belt (50–55°S; Fig. 1), making it ideally suited to study past and current changes in the SHW, temperature, and precipitation. It has a harsh, cool, wet, oceanic climate with low seasonality and high wind velocities throughout the year. Together, these represent the influence the SHW have on climate in the region (Selkirk et al., 1990). The SHW prevail almost exclusively from the west and north-west with a mean annual wind speed of 35 km h⁻¹ and gusts reaching 185 km h⁻¹ (between 1948–2025; BOM 2025). The continual dominance of the SHW drive environmental and ecosystems responses on a west-to-east gradient across the island (Chau et al., 2019; Meredith et al., 2022), including the deposition and accumulation of wind-blown inputs such as sea spray and minerogenic aerosols (Buckney and Tyler, 1974; Saunders et al., 2009). Mean annual temperature ranges from 3.1–6.6 °C, and the island experiences high annual rainfall of > 1000 mm with > 317 rainy days yr⁻¹ (between 1948–2025; BOM 2025). Rainfall has increased in recent decades with a higher frequency of intense rainfall events, mostly occurring during winter, which are then accompanied by drier, windier summers (Andersen et al., 2009; Kong et al., 2025). Persistent cloud cover over the island results in low light levels and sunshine hours per day (average 2.4 h from 1948–2022; BOM, 2025).

Figure 1 Location of Macquarie Island in the Southern Ocean and mean annual wind speeds around the Southern Hemisphere (ERA5 reanalysis data 1960–2025), showing that the island lies within the modern core Southern Hemisphere westerly wind belt (50–55° S).

Macquarie Island is geologically unique, being the only location worldwide where an intact marine ophiolite sequence of oceanic crust and upper mantle is exposed above sea level (Davis, 1987). The island is composed mostly of pillow basalts with interspersed flows of massive basalt (Selkirk et al., 1990). Dolerite, ultrabasics and intrusives are also present but are confined to the northern third of the island (Mawson, 1943). As widespread glaciation did not occur during the Last Glacial Maximum (26–20 ka), marine, periglacial and subaerial erosional, rather than glacial processes, shaped the island as well as lake formation and ontogeny. The island is fringed by a low coastal terrace leading to steep-sided slopes (20–40°) that rise to form the island plateau sitting at ∼ 200–400 m a.s.l. (Selkirk et al., 1990; McBride and Selkirk, 1998).

Figure 2 Maps showing lake sites across Macquarie Island, (a) Diatom surface sediment sites (coastal sites were originally sampled by Saunders et al. (2009); (b) Lake water chemistry sites. Colours show lake types, based on dominate hydrogeochemical processes, identified by Meredith et al. (2022). Black dots indicate that the site is also included in the diatom dataset.

The island has numerous shallow and deep lakes and ponds across the plateau and coastal terrace (Fig. 2). High accumulation of surface water and a high water-table at or very near the surface lead to the formation of extensive mires across the island (Löffler, 1984). While lake edges can form thick ice cover during winter, complete freezing of the lakes is not typically observed (Evans, 1970; Selkirk-Bell and Selkirk, 2013). The island is vegetated by bryophytes, tussock grass, herbs and sedges, with no shrub or tree species present (Selkirk et al., 1990).

2.2 Data collection

Surface sediments and water samples were collected from lakes on Macquarie Island during the 2022–2023 austral summer (referred to as 2022). Sites were selected to replicate those sampled in 2018 that were published by Meredith et al. (2022).

Lake surface sediments were collected from 30 plateau (inland) sites for diatom analyses (lake ID = LK), representing conditions more than 10 years post-rabbit and rodent eradication. Surface sediments (top 2 cm) were collected from each site using a long-handled scoop from < 1.5 m water depth. This method for sediment collection was selected for its logistical feasibility. Sediment mixing was minimised by visually assessing sampling depth and subsampling where necessary to retain only the upper ∼ 2 cm. Sediments were generally well consolidated and remained intact during collection. Based on available lead-210 (²¹⁰Pb) chronologies from Macquarie Island lake cores, this interval represents approximately 10 years of accumulation (Saunders et al., 2013; Saunders et al., 2018), comparable to surface sediment sampling approaches used in previous studies (e.g. Saunders et al., 2009). An additional 17 coastal and five plateau sites sampled in 2006 (lake ID = S; Saunders et al., 2009) were included in the diatom dataset to extend the EC and nutrient gradients of the updated dataset, totalling 52 samples (Fig. 2a). Two lakes were replicated in the 2006 and 2022 seasons (S9 = LK40, S18 = LK2; Fig. 2a).

Lake water general parameters were measured in-situ at each site, including temperature, EC, dissolved oxygen (DO), and pH using a YSI ProQuatro Multiparameter Meter, with calibration performed prior to every sampling trip. DO was calibrated in water-saturated air following YSI manufacturer protocols. EC was calibrated using a 1413 µS cm⁻¹ standard solution, and pH was calibrated using a three-point procedure with pH 4.0, 7.0, and 10.0 buffer solutions. Water samples were collected at all diatom sampling sites to measure total oxidised nitrogen (TON), phosphate (${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$), and silica (Si). Additional water samples were collected from 40 plateau sites for water chemistry analysis (Fig. 2b), including major ions and stable water isotopes (oxygen [δ¹⁸O] and hydrogen [δ²H]). Each site was sampled three times across the 2022–2023 season (November, December to January, and February). All water samples were collected from ∼ 20 cm below the water surface and were filtered in-situ with 0.45 µm polyethersulpone filters into High Density Poly-Ethylene (HDPE) bottles following the method described by Meredith et al. (2009). Water samples were refrigerated (4 °C) until analysis.

Major ions, and oxygen (δ¹⁸O) and hydrogen (δ²H) stable isotopes were analysed at the Australian Nuclear Science and Technology Organisation (ANSTO). Cations and anions were analysed using inductively coupled plasma-atomic emission spectrometry (ICP-AES). δ¹⁸O and δ²H stable isotopes were analysed with a Picarro L2130-i Cavity Ring-Down Spectrometer. Values were reported as per mill (‰) deviations relative to the international standard V-SMOW (Vienna Standard Mean Ocean Water), with a reproducible precision of ± 0.2 ‰ and ± 1.0 ‰, respectively.

Nutrient data from 2006 samples (filtered at 0.45 µm) measured soluble reactive phosphate (SRP), TON (nitrate [NO₃] + nitrite [NO₂]), and silicate (Si) using an Alpkem Autoanalyser (Continuous Flow Solution Analyser), representing the operationally defined dissolved inorganic (and therefore readily bioavailable) fractions of total N, P, and Si. In contrast, the 2022 dataset measured TON, ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, and Si ions using ICP-AES at ANSTO on filtered, undigested waters, with results reported as the corresponding inorganic species and concentrations consistently at or near detection limits. Despite methodological differences, the two approaches yield consistently low and broadly comparable concentrations in the two replicate lakes across sampling trips: for S9/LK40, ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ concentrations were 0.002 mg L⁻¹ (autoanalyser) and < 0.01 mg L⁻¹ (ICP-AES), and TON concentrations were 0.006 mg L⁻¹ (autoanalyser) and < 0.06 mg L⁻¹ (ICP-AES); for S18/LK2, ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ concentrations were 0.004 mg L⁻¹ (autoanalyser) and < 0.01 mg L⁻¹ (ICP-AES), and TON concentrations were 0.005 mg L⁻¹ (autoanalyser) and < 0.06 mg L⁻¹ (ICP-AES). Furthermore, if the 2022 analyses targeted the same reactive fractions measured in 2006, the results would still fall below or close to detection limits.

2.3 Diatom preparation and identification

Diatom preparation followed methods described by McBride (2009). Cleaned diatom solutions were mounted onto slides using Norland Optical Adhesive 61. At least 300–400 frustules were counted per sample, using Differential Interference Contrast (DIC) and oil immersion at 1000× magnification on a Zeiss Axioskope microscope, mounted with a TOUPTEK camera (U3CMOS). Species identification was primarily based on sub-Antarctic taxonomy described in Van de Vijer et al. (2002); Sterken et al. (2015); Sabbe et al. (2019); and Van de Vijver (2019). Species were photographed and documented (see Supplement for an illustrated species catalogue).

2.4 Statistical analyses

2.4.1 Water chemistry

The lake water chemistry dataset was comprised of in-situ general parameters, major ion concentrations and δ¹⁸O and δ²H values from 2018 (Meredith et al., 2022) and 2022 (this study) to understand temporal variation across the island. Data from 2018 (January to February) are referred to as sampling event 1 (E1) and sampling from the 2022 season as E2 (November), E3 (December to January), and E4 (February). Shapiro–Wilk tests showed that isotope data were normally distributed (p > 0.05), whereas general parameters and ion concentrations deviated significantly from normality (p < 0.05). Consequently, parametric tests (ANOVA, t-test, Tukey's HSD) were applied to normally distributed isotope data, and non-parametric tests (Kruskal–Wallis, pairwise Wilcoxon) to non-normally distributed general parameters ion data. Principal Component Analysis (PCA) with z-score standardised data was performed to explore relationships between variables and assess the consistency of lake types identified by Meredith et al. (2022; Fig. 2b).

2.4.2 Diatom model

Diatom inference models were developed using diatom and environmental data collected from 2006 and 2022. Ordination methods were used to describe variation in the diatom dataset, explore diatom-environment relationships, and identify unique variance explained by environmental variables. Environmental variables included were EC, temperature, DO, pH, TON, ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, and Si, with mean 2022 values used. Additional major ions were not included as these data were not available for 2006 sites. Water depth was not included as a variable as all sediment samples were collected within a narrow range of 0–1.5 m water depth, rendering ecological changes in depth negligible. Additionally, water depth is often regarded as a composite variable that acts as a surrogate for complex environmental gradients (e.g., habitat type, light, salinity, nutrients, oxygen, and taphonomy) that are largely unknown and unquantified, and therefore its inclusion can lead to spurious and misleading results (Birks et al., 1998; Juggins, 2013). Weighted Averaging (WA) was applied to determine species ecological optima and tolerance. Together WA, Weighted Averaging Partial Least Squares (WAPLS), and Maximum Likelihood (ML) models were used to develop diatom transfer functions, with cross-validation used to assess model robustness.

The relative abundance of each diatom species in each sample was calculated as the percentage of the total number of frustules counted per sample. Species occurring at ≤ 1 % relative abundance were excluded from the dataset. A full species list can be found in the Supplement. Nutrient values that were below the limit of detection were substituted with the respective detection limit value (${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, = 0.01 mg L⁻¹, Si = 0.1 mg L⁻¹, TON = 0.06 mg L⁻¹). Environmental variables were screened for skewness, with temperature, EC, ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, Si, and TON log (x+1) transformed.

PCA was performed on transformed environmental data to identify the primary gradients of environmental variation across sites. Detrended Correspondence Analysis (DCA) with detrending by segments and downweighting of rare species was performed on untransformed species data to determine whether species distributions were linear or unimodal. As the DCA axis 1 gradient length (8.2 deviation units) was > 4, unimodal ordination methods were deemed appropriate (Ter Braak and Prentice, 1988). Species data were log (x+1) transformed for remaining analysis.

A series of Canonical Correspondence Analyses (CCA) were then performed with forward selection, and scaling focused on inter-species distances, biplot scaling and downweighting of rare species. Variance Inflation Factors (VIF) of environmental variables were used to assess collinearity. As no variables had a VIF > 10, none were excluded. A full CCA, with all environmental variables included, was first performed to quantify the total amount of species–environment variance explained by the full set of variables. A series of independent and partial CCAs with variance partitioning were performed to constrain analyses, assess the relative explanatory power, and assess the unique and shared variance contributions of each variable. Individual CCAs, of each variable alone, estimate the marginal (unconstrained) explanatory power (i.e., how much variation a single variable explains when considered alone, without accounting for correlations with other variables). Partial CCAs assess the unique (conditional) contribution of each environmental variable after statistically controlling for all remaining variables. This analysis isolates the variance uniquely attributable to each predictor and identifies variables whose explanatory power is driven by covariation with others. Finally, variance partitioning was used to decompose the total explained variation into unique and shared fractions, allowing assessment of how much variation was due to individual predictors versus overlapping environmental gradients. Permutation test results (p > 0.05), CCA coefficients and lambda ratios (${\mathit{\lambda }}_{\mathrm{1}/}{\mathit{\lambda }}_{\mathrm{2}}$) of the first constrained eigenvalue (λ₁) to the second unconstrained eigenvalue (λ₂) were used to identify the environmental variables most appropriate for quantitative inference models. As a guide, high λ₁/λ₂ ratios are necessary for a variable to have enough explanatory power to be included in quantitative inference models (Ter Braak and Prentice, 1988; Juggins, 2013). All ordination analyses were performed using the vegan package version 2.7–1 (Oksanen et al., 2013) in R (R Core Team, 2024).

ML and iterations of inverse (_INV) and classical (_CLA) WA models with and without tolerance downweighing, and WAPLS with up to five components were assessed to determine the best performing transfer functions. These methods were applied because they capture different aspects of species–environment relationships: WA provides a simple unimodal estimator; WAPLS allows more complex responses through latent components; and ML emphasises taxa with narrow ecological tolerances. Using multiple approaches therefore offers complementary strengths and helps identify the most reliable and robust model through cross-validation. All models were performed with bootstrapping and 100 iterations. Model R², bootstrapped R² ($R_{\text{boot}}^{\mathrm{2}}$), root mean square error (RMSE) and root mean square error of prediction (RMSEP) values were used to assess performance. RMSEP and $R_{\text{boot}}^{\mathrm{2}}$ performance were favoured over R² and RMSE. RMSEP between WAPLS components were also used to assess overfitting. WA and WAPLS-1 results are often similar as WAPLS is built upon on the same weighted-averaging framework as WA (ter Braak and Juggins, 1993). When this was the case and WAPLS components did not improve performance, WA was favoured as the most parsimonious model. Software program C2 version 1.8 (Juggins, 2003) was used to develop all transfer functions.

3 Results

3.1 Lake water chemistry

Analyses of 40 plateau lakes on Macquarie Island showed that lake water general parameters (EC, pH and DO), and nutrients, did not vary significantly (p > 0.01) across the 2022 sampling events (E2–4; Table 1). Temperature varied, being significantly lower in E2 compared to E3 and E4 (p < 0.01). Lakes were moderately acidic (pH 5.7) to alkaline (pH 9.14). Mean EC ranged from 126–261 µS cm⁻¹, with a decrease in EC from west to east across the island. Lakes were oxic (DO = 8.64–12.61 mg L⁻¹) and oligotrophic, with ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ and TON concentrations under or close to detection limits (< 0.01–0.02 mg L⁻¹ and < 0.06–0.1 mg L⁻¹, respectively). Similarly, comparison with data from plateau lakes sampled in 2018 showed no significant difference (p > 0.01) across all lake water general parameters, excluding temperature, indicating generally stable conditions in plateau lakes across years. However, a comparison between plateau lakes measured in 2022 and coastal lakes in 2006 did show significant differences (p < 0.01). Coastal sites in 2006 were generally eutrophic with higher nutrient ranges (TON = 0.007–4.636 mg L⁻¹ and ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ = 0.1–9.9 mg L⁻¹), and higher EC (406–1482 µS cm⁻¹), while temperature, DO, pH, Si were not significantly different (p > 0.01) (Table 1; see Table S1 and S2 in the Supplement for full results).

Table 1 Summary of lake water general parameters and nutrient data. Temp. = temperature, DO = dissolved oxygen, EC = electrical conductivity, Si = silicate, ${\mathrm{PO}}_{\mathrm{4}}^{-\mathrm{3}}$ = phosphate, TON = total oxidised nitrogen.

Figure 3 Box and whisker plots showing the range and mean of 2018 and 2022 lake water major ions (a) Cl and Na; (b) SO₄, Br, Ca, Mg, K, SiO₂, Fe, F, and Al; and (c) stable water isotopes δ²H and δ¹⁸O.

Major ion analysis of the 40 plateau lakes showed that, although dilute in concentration, Cl (1.1–3.7 mmol L⁻¹) and Na (0.9–2.6 mmol L⁻¹) dominate the ionic composition of all lake waters (Fig. 3; see Table S3 for full cation and anion results). All lakes showed similar ionic ratios to seawater for SO₄, Cl, Mg, and Na, suggesting a marine origin. Seawater ionic ratios diverged for K, Ca and F for some lakes, while SiO₂ was higher in all lakes, suggesting additional sources for these ions (Fig. 4).

Figure 4 Scholler plot comparing the ionic composition (SO₄, Cl, Mg, Na, K, Ca, F and SiO₂) of Macquarie Island lake waters and seawater (red line), categorised into sampling years 2018 (black lines) and 2022 (grey lines).

Statistical analyses showed no significant (p > 0.005) differences in major ion concentrations across 2022 sampling events (E2–4). Broader changes were detected between 2018 and 2022, with Cl, SO₄, Br, and Mg all showing significantly higher mean concentrations (p < 0.005) in 2018 compared to all 2022 sampling events (E2–4). Fe, Na, K, Ca, F, SiO₂, and Al did not significantly vary (p > 0.005) between sampling events. All ions that show significant variation (p < 0.005) have predominantly marine sources.

Figure 5 Principal Component Analysis (PCA) of Macquarie Island lake waters, showing the relationships between major ions and environmental parameters. Lakes are coloured by lake type, showing that lakes cluster based on the dominate geochemical processes identified by Meredith et al. (2022). Dist. west coast = distance from the west coast (m), Conduct = electrical conductivity.

A PCA shows the relationship between lake water chemistry parameters, with samples grouped by lake type and sampling year (Fig. 5). Together, PC1 and PC2 captured 53 % of the total variance in the dataset. PC1 represents a salinity and sea-spray gradient with variability in EC, distance from the west coast, Na, Cl, Br, Ca, Mg, and K captured. PC2 represents an altitude and terrestrial ion gradient with variability in elevation, temperature, SiO₂, Ca, Fe, and F captured. Lakes cluster according to environmental processes (groups derived from Meredith et al., 2022), with SSA influenced lakes having positive PC1 scores, which suggests higher concentrations of marine derived ions. The grouping of samples influenced by catchment processes and those influenced by rainfall is driven by PC2 with catchment influenced lakes having lower elevation and higher ion concentrations. SSA and rainfall influenced lakes cluster based on the year that they were sampled with greater ion concentrations in 2018.

3.2 Stable water isotopes

δ²H and δ¹⁸O were measured in the 2022 samples and ranged from −38.3 ‰ to −9.2 ‰ and −38.3 ‰ to −0.67 ‰, respectively. Lake waters in 2018 and 2022 fell below the Global and Cape Grim (northwest Tasmania) Meteoric Water Lines (MWLs), suggesting slight isotopic enrichment in Macquarie Island's lakes (Fig. 6). In 2018, lake waters were significantly higher in δ²H (mean −24.8 ‰) and δ¹⁸O (mean −2.83 ‰) compared to 2022 (mean δ²H = −24.8 ‰ and δ¹⁸O = −3.55 ‰). Significant isotopic enrichment of δ²H and δ¹⁸O (p < 0.001) can be seen in the data at the beginning of the 2022 austral summer (E2–3; Fig. 6a). SSA influenced lakes tended to have higher isotopic values, while catchment influenced lakes had lower values (Fig. 6b), with a significant difference between all lake types detected (p < 0.001). LK20 and LK21 from E2 were outliers, plotting above the MWLs with lower δ¹⁸O values. Cl concentrations appeared to be related to δ²H and δ¹⁸O values (Fig. 6c), however the correlation between the parameters was low (R² ≤ 0.24). This lack of relationship was consistent across lake types and sampling events.

Figure 6 Stable water isotope δ²H and δ¹⁸O differences in Macquarie Island lake waters, shown across: (a) sampling events; (b) lake type; (c) δ²H and δ¹⁸O relationship shown with Cl concentration (mmol L⁻¹). Solid black line shows Macquarie Island regression, blue dashed line is the Global Meteoric Water Line (GMWL: δ²H = 8δ¹⁸O+10), and black dashed line is the Cape Grim Local Meteoric Water Line (LMWL: δ²H = 6.8δ¹⁸O+6.65).

3.3 Diatoms

3.3.1 Diatom communities

In total, 141 diatom taxa from 45 genera were identified in the 52 plateau and coastal lakes. 96 taxa (including 21 unknown species) from 30 genera remained in the dataset after taxa with ≤ 1 % relative abundance were excluded. The most taxon-rich genera were Pinnularia (16 taxa), Psammothidium (12 taxa), and Planothidum (8 taxa). The most dominant taxa being both common, occurring in > 15 lakes, and abundant, occurring > 25 % relative abundance in at least one sample, were Aulacoseira principissa Van de Vijver, Psammothidium abundans (Manguin) Bukhtiyarova and Round, Psammothidium confusum var. atomoides (Manguin) van de Vijver, unknown species 111, Psammothidium confusum (Manguin) van de Vijver, Fragilaria capucina Desmazières, and Navicula bergstromiana Vyverman et al. Coastal and plateau lakes showed distinctly different assemblages, with coastal lakes exhibiting less diversity (mean number of species = 20) and dominated by taxa including F. capucina, unknown sp. 111, Planothidium quadripunctatum (D. R. Oppenheim) Sabbe, Planothidium delicatum (Kützing) Round and Bukhtiyarova, and Planothidium lanceolatum (Brébisson ex Kützing) Lange-Bertalot. Plateau lakes were more diverse (mean number of taxa = 40) and dominated by A. principissa, P. abundans, P. confusum var. atomoides, P. confusum, N. bergstromiana, Achnanthidium modestiformis (Lange-Bertalot) Van de Vijver, Cocconeis placentulata Ehrenberg, and unknown species 21. No taxa were found in all lakes, and none were uniquely restricted to either coastal or plateau lakes, although clear differences in species composition and relative abundance were observed. See Fig. 7 for microscopy photos of the most abundant taxa.

Figure 7 The most abundant diatom taxa from plateau and coastal lakes on Macquarie Island. (a) Aulacoseira principissa; (b) Psammothidium abundans; (c) Psammothidium confusum; (d, e) Cocconeis placentulata; (f) Navicula bergstromiana; (g) Fragilaria capucina; (h) Unknown species 21; (i) Psammothidium confusum var. atomoides; (j) Achnanthidium modestiformis; (k) Unknown species 111; (l) Planothidium delicatulum; (m) Planothidium lanceolatum.

3.3.2 Diatom-environment relationships

The full CCA model was significant (p = 0.001), with environmental variables explaining 23.9 % (Table 2) of the total variance in diatom species composition (constrained inertia = 1.5). Together, the first two canonical axes explained 62.8 % of the total constrained variance (Table 2; Fig. 8). Forward selection with Bonferroni corrections, identified EC, pH (p = 0.003), and temperature (p = 0.015) as the most significant predictors of diatom community composition, collectively explaining 15.4 % of the total variance. This equates to 70 % of the total explained variance in the full model, capturing the major environmental gradients influencing species distribution with a more parsimonious model.

Table 2 Full CCA model results.

Figure 8 Full CCA ordination biplot of diatom species and environmental data, numbers indicate sites, and red symbols indicate diatom species. Si = silicate, P = phosphate, N = total oxidised nitrogen, and conduct = electrical conductivity (EC).

Table 3 Individual CCA results, independent CCAs run for each variable. (^∗ = significant p value < 0.05).

Individual CCAs were performed to assess the total explanatory power of each variable. EC and pH were shown to be the strongest, individually explaining 10.62 % and 5.43 % of the total variation, respectively. This corresponds to 44.4 % and 22.7 % of the total variance in the full CCA model. Additionally, EC was the only variable with a high $\mathit{\lambda }\mathrm{1}/\mathit{\lambda }\mathrm{2}$ ratio (${\mathit{\lambda }}_{\mathrm{1}}/{\mathit{\lambda }}_{\mathrm{2}}$ = 1.31; Table 3), suggesting it is the only variable with enough explanatory power for inference modelling.

Table 4 Partial CCA results, where each variable was tested with the covariation of other variables controlled (^∗ = significant p value < 0.05).

Partial CCAs were performed with each environmental variable tested separately while controlling for covariation with all other variables, to quantify unique and shared variance contributions. EC, pH, and Si were the only variables to have significant unique contributions (p ≤ 0.01; Table 4). The shared and unique variance of each environmental variable is shown in Fig. 9. EC explained the largest proportion of total constrained variation (46 %) in diatom community composition, with a large shared component (18 % unique, 28 % shared), suggesting it acts along a major environmental gradient shared with TON and ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ (Fig. 8). Despite this, it performed well in all other CCAs and its unique contribution remained high, indicating it is an important independent driver of diatom structure across Macquarie Island lakes. Furthermore, low VIFs among all environmental variables (VIFs < 3) indicated that multicollinearity was low. EC (VIF = 2.6) showed a low correlation with other variables (R² ≤ 0.47), suggesting it represents a largely independent gradient in the dataset. In contrast, pH had similar unique variance (18.4 %) and lower shared variance (2.4 %), implying a more independent ecological influence.

Figure 9 Variance partitioning showing unique and shared proportions of variance explained by each environmental variable.

3.3.3 Species optima and tolerances

Species optima across major environmental gradients EC and pH were determined with WA. F. capucina, P. lanceolatum, and P. delicatulum were found across the EC range with broad tolerances, however each species showed different optima (Fig. 10). The apparent bimodal distribution of F. capucina likely reflects ecological plasticity across differing hydrochemical conditions and/or potential taxonomic aggregation within this morphotype. Exploration of a GAM-based response curve (Fig. S2 in the Supplement) indicates F. capucina has a weak non-linear relationship with EC, suggesting a broad and flexible ecological response rather than a strongly defined unimodal optimum. Unknown sp. 111 was found to tolerate high EC, while most other dominant species, including A. principissa, C. placentulata, N. bergstromiana and dominant Psammothidium species show tolerance and optima for low EC.

Figure 10 Weighted averaging (WA) optima (solid red line) and tolerance ranges (dashed red lines) for dominant diatom species across the electrical conductivity gradient in Macquarie Island lakes. Observed relative abundances (%) are plotted against log-transformed electrical conductivity, with fitted loess curves illustrating species response shapes along the conductivity gradient (126–1482 µS cm⁻¹). Roman numerals indicate Gaussian response curve type.

Figure 11 Weighted averaging (WA) optima (solid red line) and tolerance ranges (dashed red lines) for dominant diatom species across the pH gradient in Macquarie Island lakes. Observed relative abundances (%) are plotted against pH, with fitted loess curves illustrating species response shapes along the pH gradient (5.50–9.14 °C).

Dominant diatom species occurred across the pH gradient, with species optima ranging from moderately acidic to neutral (Fig. 11). Most dominant species had optima for moderate acidity, with Psammothidium species, A. principissa, and N. bergstormiana showing narrow pH tolerances. C. placentulata and F. capucina showed broader tolerances and unknown sp. 111 was the only dominant species with a near neutral optima.

3.3.4 Diatom transfer functions

Ordination analyses showed that EC and pH explained significant and independent proportions of variance in diatom composition. While temperature also showed significant but lesser contributions, it was not considered for transfer function development as diatom-based temperature reconstructions as species responses to temperature can be indirect and influenced by multiple co-varying environmental gradients, limiting the reliability of temperature inference (e.g. Juggins, 2013). While Si was shown to independently contribute to diatom variance, reduced CCA modelling with forward selection did not indicate it to be a major environmental gradient. Transfer functions were therefore only developed for EC and pH. Transfer function results for the best performing WA, WAPLS, or ML model for EC and pH are described in Table 5.

Table 5 Best performing WA, WAPLS or ML model results for electrical conductivity and pH.

For EC, WA_INV and WAPLS-1 produced near identical results. WA_INV was favoured as the simpler model (WA_INV, R² = 0.83, $R_{\text{boot}}^{\mathrm{2}}$ = 0.74, RMSE = 0.31, RMSEP = 0.39). Although WAPLS-2 to -5 increased R² and reduced RMSE, each successive component progressively increased RMSEP by 13 %–16 %, thereby reducing performance. Given the unimodal gradient structure of the dataset, ML modelling was also assessed for EC. ML showed slightly stronger predictive performance to WA_INV, with higher $R_{\text{boot}}^{\mathrm{2}}$ and comparable RMSEP (ML, R² = 0.91$,R_{\text{boot}}^{\mathrm{2}}$ = 0.80, RMSE = 0.23, RMSEP = 0.40).

Figure 12 Comparison of observed environmental measurements with values predicted by diatom-based transfer functions: (a) pH estimated using WAPLS-2; (b) conductivity estimated using WA with inverse deshrinking (WA_INV); and (c) conductivity estimated using the maximum likelihood (ML) method. Black lines show the 1:1 line.

Comparison of observed and predicted value scatter plots indicated that ML achieved a tighter fit, with WA_INV showing increased predictive error at higher EC ranges (Fig. 12). However, further inspection indicated that only ∼ 30 % of taxa displayed Gaussian (Types IV–V) response curves indicating ML may not be the most appropriate approach (Fig. 10 and see Table S5 for full Gaussian response curve results). However, overfitting from ML is not likely as RMSEP did not increase. Overall, both models show cross-validated performance and are considered robust.

pH had poor performance across all WA and WAPLS models with high RMSEP (≥ 0.6) and low $R_{\text{boot}}^{\mathrm{2}}$ (≤ 0.26). WAPLS-2 was found to be the strongest model (WAPL-2, R² = 0.72, $R_{\text{boot}}^{\mathrm{2}}$ = 0.26, RMSE = 0.33, RMSEP = 0.60). ML modelling showed poor performance for pH with high RMSEP = 0.93.

3.4 Discussion

3.5 Annual and seasonal lake water hydrogeochemical variation

The seasonal hydrochemistry dataset presented in this study support the 2018 baseline assessment (Meredith et al., 2022) that lake water chemistry is controlled by SSAs, terrestrial catchment processes, elevation and rainfall dilution. The seasonal water chemistry data, which is presented for the first time in this study, shows that there is no significant variation in major ions across the 2022–2023 austral summer. Comparison of 2018 and 2022 data shows that significant variation (p < 0.005), in some major ions is evident, with higher concentrations in 2018 of Br and Cl, SO₄, and Mg associated with SSAs. Although not statistically significant (p > 0.005), other sea-spray derived ion mean values (Na, Ca, K) were higher in 2018 (Fig. 3; Table S4). However, not all ions increased in concentration, and terrestrially derived ions such as Fe, F, Si, and Al were lower in 2018, suggesting that lakes in 2022 had stronger SSA influence.

Despite these changes, PCA of the 2018 and 2022 lake water chemistry datasets (Fig. 5) show that lakes typically cluster by the lake types identified in Meredith et al. (2022) (i.e., as SSA, catchment and rainfall influenced). This indicates that the major hydrogeochemical processes influencing Macquarie Island lakes are consistent between years, with no major environmental shifts occurring in weathering and erosion of the island's geology, suggesting these processes are stable on an annual time scale. Identifying hydrogeochemical stability is important for identifying lake sites suitable for diatom–conductivity inference models. It also strengthens palaeoclimate interpretations by suggesting that local lake dynamics are relatively constant, supporting the hypothesis that in sea-spray dominated lakes, proxies primarily record externally forced changes driven by the SHW rather than internal hydrological or geochemical dynamics (Saunders et al., 2018; Perren et al., 2020). This further emphasises that water chemistry characteristics are critical to consider in site-selection to develop reliable SHW reconstructions on Macquarie Island.

3.6 Evaporation

Interpreting environmental proxies as direct indicators of climate variability can be challenging, as multiple processes may produce similar signals (Molén, 2024). When using diatom–conductivity models to infer past SHW variability, it is essential to consider how near surface evaporation has the potential to concentrate ions in surface waters and mimic the effects of other processes such as increasing lake water salinity from SSA deposition due to stronger winds. Although this study and previous studies (Evans, 1970; Buckney and Tyler, 1974; Meredith et al., 2022) have demonstrated that SHW-driven SSA inputs are a dominant control on lake water chemistry on Macquarie Island, the role of evaporation in amplifying these signals remains unclear.

To explore this, we analysed δ²H and δ¹⁸O values from 2018 and across the 2022–2023 Austral summer. Isotopic enrichment is evident across the 2022–2023 season (Fig. 6a), indicating a strengthening evaporative signal through summer. Lake stable water isotopes sampled in 2018 were significantly more enriched than the 2022 mean, as expected given that the 2018 samples were collected in late summer, when evaporative effects are strongest and cumulative due to warmer temperatures throughout summer. This is supported by lake water temperatures being significantly lower in early summer (E2) compared to late summer (E3 and E4; Table 1). Given Macquarie Island's persistently high cloud cover, humidity, and low sunshine hours (BOM, 2025), solar evaporation is likely limited and confined to the summer season. Evaporation can produce heavy-isotope enrichment and the residual lake water becomes progressively enriched in heaver isotopes (δ²H and δ¹⁸O), moving away from the Global MWL (Gat, 1996). Comparisons between E1 (2018) and E4 (2022), which were sampled at the same time of the year in January–February provide a valuable comparison of potential interannual variability in lake water chemistry and processes (Fig. 3, Table S4). These two sampling events show near identical mean isotopic composition (p > 0.01; δ²H = −20.7 ‰ and δ¹⁸O = −2.8 ‰ in 2022, and δ²H = −20.1 ‰ and δ¹⁸O = −2.8 ‰ in 2018), suggesting broadly stable summer evaporative conditions between years. Furthermore, SSA influenced lakes on the plateau are in proximity to the west coast and have the greatest exposure to the SHW (Fig. 2b). These lakes have significantly higher isotopic enrichment (Fig. 5b), providing further evidence that wind is likely the primary driver of evaporation in plateau lake waters across Macquarie Island, particularly in lakes located on the west coast, which may be most suitable for reconstructions of SHW dynamics (Saunders et al., 2018). As both wind-enhanced evaporation and wind-driven SSA transport and deposition contributes to the concentration of SSA ions in these lakes, both ion deposition and concentration reflect a SHW signal.

Cl⁻ is a robust tracer of hydrogeochemical processes and, together with δ²H and δ¹⁸O values, can be used to better understand evaporation (Kirchner et al., 2010). While the isotopic enrichment observed generally indicates an evaporative signal, the absence of a correlation between Cl⁻ and δ²H or δ¹⁸O (Fig. 6c) suggests that isotopes are capturing short-term (summer) evaporation rather than sustained evaporative concentration sufficient to increase Cl⁻ concentration in the lake waters like those in environments driven primarily by solar evaporation (Meredith et al., 2009). On Macquarie Island, wind-driven SSA deposition and rainfall dilution therefore likely remain the primary drivers of Cl⁻ variability. Consequently, EC in lake waters of Macquarie Island remains a robust proxy for interpreting variations in SHW strength.

3.7 Diatom communities

Diatom analysis showed that typical sub-Antarctic genera (Van de Vijver, 2019; Goeyers et al., 2022), including Psammothidum, Planothidium, and Fragilaria, dominated lake diatom communities on Macquarie Island. Across 52 lakes, 141 taxa were identified, indicating intermediate species diversity relative to previous studies on the island which reported 102 (McBride, 2009) and 208 (Saunders et al., 2009) species. Consistent with these earlier studies we have demonstrated that diatoms on Macquarie Island exhibit clear and distinct ecological preferences. Combined with the pronounced environmental gradients among lakes, these species-environment relationships provide a strong basis for using diatoms as indicators of limnological conditions and environmental change.

Psammothidium species are characteristic of low EC sites (Van de Vijer et al., 2002), and while abundance and dominance of key Psammothidium species showed variation along the lower end of the EC gradient (Fig. 10), they dominated low EC sites. N. bergstromiana is considered endemic to Macquarie Island and was commonly found with dominant Psammothidium species at low EC sites, typically occurring where EC was < 200 µS cm⁻¹, consistent with what has previously been reported (Sabbe et al., 2019). A. principissa, previously identified as Aulacoseria distans (Ehrenberg) Simonsen on Macquarie Island (McBride, 2009; Saunders et al., 2009), was also common at low EC (Fig. 10). This taxa is commonly found on sub-Antarctic Islands and is suggested to prefer very low conductance values < 80 µS cm⁻¹ (Van de Vijver, 2012), while EC was not observed < 160 µS cm⁻¹ in this dataset, A. principissa may be an indicator of very low EC conditions.

P. lanceolatum was a dominant high EC, high nutrient taxa. While it has been found to dominate flora elsewhere, this contrasts previous studies where it has been reported to be characteristic of oligotrophic conditions (Van de Vijer et al., 2002). F. capucina, P. delicatulum, and unknown sp. 111 were more commonly dominant at high EC. F. capucina was found across the EC gradient displaying a bimodal distribution (Fig. 10), consistent with its cosmopolitan and ecologically tolerant nature (Van de Vijer et al., 2002).

3.8 Developing transfer functions

A key aim of this study was to update and improve existing quantitative diatom models for Macquarie Island. While the dominant taxa identified here are consistent with those reported by Saunders et al. (2009), the strength of some diatom–environment relationships differ. EC and pH remain strong explanatory variables for diatom variation, whereas ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$ and Si showed limited influence in the present study (Tables 3 and 4). This reduced explanatory power likely reflects the low nutrient variability across plateau lakes, where concentrations were generally below detection limits.

Sub-Antarctic lakes are characteristically oligotrophic; high nutrient levels do occur, but they are associated with peatlands or animal colonies, as is the case with coastal lakes on Macquarie Island (Selkirk et al., 1990). Saunders et al. (2009) recorded greater nutrient variability across plateau sites at the lower end of the EC gradient, attributed to enhanced organic inputs during periods of high ecological disturbance from invasive rabbits. These differences suggest that the dataset in the present study represents post-eradication limnological conditions that are more reflective of pre-invasion or near-natural states, providing an updated basis for developing robust diatom–environment models that are less influenced by disturbance. The use of field parameters collected across multiple sampling events increases confidence that the developed models reflect representative environmental conditions (Goldenberg Vilar et al., 2018; Kennedy and Buckley, 2021). Such repeated-sampling approaches are rarely achieved in diatom transfer function development, particularly in remote regions such as the sub-Antarctic.

Widespread recovery of vegetation communities provides evidence for catchment scale ecosystem recovery across the island (Springer, 2018; Fitzgerald et al., 2021). While quantitative runoff or nutrient time series are not available to directly extend this to limnological conditions, we are able to provide site-specific evidence from one sedimentary diatom record. At Emerald Lake (LK6), downcore diatom assemblages show a clear ecological shift coincident with the introduction of rabbits (1878 CE) and their establishment on Macquarie Island, with F. capucina and P. abundans dominating downcore intervals (Saunders et al., 2013). In contrast, diatom assemblages in recent (2022) surface sediments from this site exhibit higher diversity (48 species) and greater similarity in assemblage composition to pre-rabbit sediment intervals rather than assemblages from previously collected 2006 surface sediments (15 species), which were dominated by F. capucina (48 % relative abundance). Notably, F. capucina was absent from the 2022 surface sample. Together, these lines of evidence suggest that the modern calibration dataset is less influenced by organic inputs and erosion associated with the rabbit invasion period. Accordingly, we interpret the 2022 dataset as representing post-eradication recovery conditions that are moving toward, but do not necessarily fully reflect, pre-disturbance baseline states. However, further studies are needed to assess how widespread this is across the island.

EC was shown to be the major independent driver of diatom assemblages, with strong performance in all CCA models, individually accounting for almost 50 % of the variance in the full CCA model and strong explanatory power as indicted by ${\mathit{\lambda }}_{\mathrm{1}}/{\mathit{\lambda }}_{\mathrm{2}}$ (Table 3). Although EC reflects a high proportion of shared variance (Fig. 9), this is consistent with its role as an integrative measure of ionic strength and catchment inputs. The shared component primarily reflects its covariation with major ions and nutrient variables (TON and ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$), which are typically correlated with EC in these systems. Despite this overlap, EC retained a strong and highly significant independent effect (p = 0.001), confirming its dominant ecological influence on diatom distributions.

While ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, TON, and Si each explained significant but moderate portions of individual variance (Table 3), they were no longer significant once covariation was controlled for (Table 4), meaning their explanatory power is mostly shared variance with other environmental gradients, primarily EC and each other, with negligible unique variance. This interpretation is supported by simple linear correlations between EC and nutrient variables in coastal lakes, which show weak or absent relationships (Fig. S3). These results indicate that although nutrients and EC co-occur in coastal systems, nutrient variability is not strongly or systematically coupled to EC. Together with the VIFs (< 3), variance partitioning and partial CCA results (Table 4), this supports the interpretation that the weak unique nutrient signal reflects an ecological reality in which EC exerts a first-order control on diatom assemblages, rather than a statistical artefact of collinearity. Similarly, while pH explained a major and independent gradient in diatom variation within plateau lakes (Fig. 8), it did not capture assemblage changes across high-EC, high-nutrient sites. This was indicated by individual CCA results, which were less than half of the variance explained by EC (Table 3). This, paired with the widespread oligotrophic nature of plateau lakes on Macquarie Island lends strength to the independent explanatory power of EC across the whole dataset.

Furthermore, pH showed poor predictive performance as a transfer function, with the lowest $R_{\text{boot}}^{\mathrm{2}}$ = 0.26, and high RMSEP = 0.6 from the WAPLS-2 model (Table 5). While this is surprising due to the strong pH gradient across plateau sites, most diatoms were found across the pH gradient (Fig. 11) with some species showing broad tolerance and pH ranges. EC had the strongest performance with the WA and ML models producing the highest $R_{\text{boot}}^{\mathrm{2}}$ (0.74 and 0.80, respectively) and comparable RMSEP. WA_INV and WAPLS-1 showed identical performance, with no benefit from additional WAPLS components, which progressively increased predictive error and decreased $R_{\text{boot}}^{\mathrm{2}}$, suggesting overfitting. WA was therefore chosen over WAPLS as the simpler model.

The WA EC transfer function performed better than the previously published Macquarie Island diatom-conductivity transfer function (Saunders et al., 2009), with higher $R_{\text{boot}}^{\mathrm{2}}$. However, some caution is warranted when predicting across the upper EC range, where greater predictive error is evident (Fig. 12). This can be attributed to lower species turnover, higher variability in ${\mathrm{PO}}_{\mathrm{4}}^{\mathrm{3}-}$, TON, and Si, and fewer sites at the upper end of the nutrient and EC gradients. Further refinement of the EC transfer functions could be achieved with more evenly distributed sampling across the environmental gradient. The ML model, with higher $R_{\text{boot}}^{\mathrm{2}}$, appears more capable of addressing these issues and maintains more consistent predictive power across the EC range. This is likely due to its explicit curve-fitting approach. By estimating individual species optima and tolerances, ML can better represent asymmetric or skewed response curves (Birks, 2012). The ML transfer function is therefore considered to be the preferred model, although both WA and ML are robust based on comparable RMSEP.

3.9 Future applications for reconstructing past climate changes

The conceptual link between large-scale wind regimes and long-term limnological and ecological responses in Southern Hemisphere lake systems is well established in previous studies (e.g., Saunders et al., 2009, 2015, 2018; Perren et al., 2020, 2025; Van Nieuwenhuyze, 2020; Humphries et al., 2021; Meredith et al., 2022). The data presented here builds on this existing framework and demonstrates that incorporating diatom data with seasonal and multi-year hydrogeochemical data provides a unique opportunity to comprehensively understand diatom-environment responses. By quantifying temporal variability in hydrogeochemical processes, including the role of evaporation, this study strengthens confidence that EC reflects SHW-driven sea-spray inputs rather than local lake hydrogeochemical processes. This hydrological context is critical for interpreting diatom–environment relationships and ensuring the reliability of EC as a proxy for past SHW behaviour, providing a strong foundation for future palaeoclimate reconstructions. The resulting diatom–conductivity model provides a robust and ecologically grounded framework for reconstructing long-term SHW variability on Macquarie Island and establishes an important benchmark for sub-Antarctic palaeoclimate comparisons across the region. This model will be applied in future studies to reconstruct past variability in the SHW and associated hydroclimatic changes on Macquarie Island. Any future time-series monitoring of lake EC paired with local wind speed records would further allow direct assessment of the relationship between wind speed and lake salinity, including the rate of hydrogeochemical response and any wind speed thresholds required to drive measurable change.

By capturing post-eradication and near-natural ecological conditions, the EC model developed in this study offers an improved foundation for assessing long-term wind-driven variability, as it reduces ecological noise associated with past disturbance. When applied in parallel with other proxies, such as isotopic or geochemical indicators, these reconstructions will contribute to a more comprehensive understanding of past SHW dynamics and their role in modulating Southern Hemisphere mid-high latitude climate, thereby providing context for understanding future changes.

Furthermore, a multiproxy approach will be valuable for independently reconstructing key climatic drivers, including precipitation, temperature, and atmospheric circulation, thereby improving interpretations of past SHW variability and helping to assess how hydroclimatic processes may modify EC signals (e.g., through dilution and enrichment). On Macquarie Island, geochemical indicators of sea-spray and dust inputs (e.g., S, Br, Ti) can help distinguish marine aerosol delivery from catchment-derived material, while glycerol dialkyl glycerol tetraethers (GDGT)-biomarker reconstructions can provide an independent constraint on temperature variability. Both approaches are currently being undertaken on Macquarie Island lake sediment cores by our research group. Mercury (Hg) concentrations and isotopes offer an independent proxies for atmospheric transport and Hg deposition linked to large-scale circulation, precipitation, and the influence of seabirds, all particularly well suited to the remote setting of Macquarie Island (Schneider et al., 2022; Guédron et al., 2018), which is also the focus of ongoing work (e.g. Schneider et al., 2024; Li et al., 2026). Although isotope (δ²H, δ¹⁸O) palaeo-records are not currently available for Macquarie Island, they represent an important avenue for future research to constrain precipitation–evaporation balance. Together, these complementary proxies provide a framework to separate the relative influence of atmospheric circulation, hydroclimate, and temperature on lake systems, providing more comprehensive palaeoclimate records and interpretations.

4 Conclusion

This study aimed to update and re-evaluate the reliability of diatom–conductivity models as a proxy for reconstructing SHW variability on Macquarie Island by analysing diatom–environment relationships in the context of seasonal and multi-year water chemistry and isotopic analyses. Our results demonstrate that although lake hydrogeochemical processes vary locally, they remain stable seasonally and between years. Lakes near the west coast and on the western edge consistently reflect strong SSA influence, and while short-term evaporative enrichment occurs during summer, it does not obscure the dominant signal of SHW-driven SSA inputs. Accordingly, EC reliably reflects SSA deposition rather than internal lake hydrogeochemical processes, providing a firm mechanistic basis for the use of EC as an indicator of SSA deposition in palaeoclimate studies on Macquarie Island.

Diatom–environment relationships were found to be strong and ecologically coherent, supporting the development of a robust diatom–conductivity transfer function. Importantly, this study highlights the need for careful site selection, with lakes that demonstrate stable hydrogeochemical behaviour, clear SSA influence, and limited local disturbance providing the most reliable archives for reconstructing past SHW variability. The resulting EC transfer function offers a reliable tool for reconstructing long-term SHW dynamics, supported by well-characterised modern hydrological controls. Together, these findings establish Macquarie Island as a well-constrained system for SHW reconstructions and provide a strong foundation for future palaeoclimate work across the sub-Antarctic region.