Hot-spring inputs and climate drive dynamic shifts in archaeal communities in Lake Magadi, Kenya Rift Valley

The methane index (MI) is an organic geochemical index that uses isoprenoid glycerol dialkyl glycerol tetraethers (GDGTs) as a proxy for methane cycling. Here, we report results from sediments in core MAG14-2A that span almost 500 ka in Lake Magadi, Kenya. The deposits show abrupt shifts between high and low MI values through calcareous, tuffaceous and zeolitic silts. The MI switches “off” (MI < 0.2) and “on” (MI > 0.5) through the core with bulk organic matter enriched in ¹³C during “MI-off” periods (∼ −18 ‰) in the upper part of the core, whereas ¹³C is lower (−22 ‰ to −25 ‰) in lower parts of the sedimentary sequence. Sediments deposited when the MI switches “on” showed δ¹³C_OM values as low as −89.4 ‰, but most values were within the range of −28 ‰ to −30 ‰, which is consistent with contributions from methanogens rather than methanotrophs. Thus, the likely source of these high MI values in Lake Magadi is methanogenic archaea. Our results show that hydrothermal inputs of bicarbonate-rich waters into Lake Magadi combined with further evaporative concentration cause a shift in the dominant archaeal communities, alternating between two stable states.

1 Introduction

Life thrives in East African soda lakes and has been the subject of modern studies of both prokaryotic and eukaryotic organisms, but few researchers have studied their sediments over geologic timescales (Schagerl, 2016, and chapters therein). Soda lakes represent ∼ 18 500 km² in East Africa (calculated from values in Melack and MacIntyre, 2016). When compared to the three largest African freshwater lakes (lakes Victoria, Tanganyika, and Malawi), these soda lakes account for ∼ 13 % of the total lake-surface area in East Africa. A survey of microbial isolate diversity in East African lakes found evidence for cyanobacterial and archaeal primary producers, with both oxygenic and anoxygenic phototrophs among the microbial population (Grant and Jones, 2016). Unique aerobic and anaerobic heterotrophs that use a variety of electron donors, including sulfur, sulfate, nitrite, carbon dioxide, and methane, were also identified (Grant and Jones, 2016, and sources therein). Moreover, many thermophilic archaea and bacteria isolates were observed near hot-spring outflows (Grant and Jones, 2016).

Saline alkaline (soda) lakes in the East African Rift often become stratified meromictic waterbodies with a dense monimolimnion below a chemocline (Melack and MacIntyre, 2016). Oxygen rarely penetrates the monimolimnion waters; as a result, anaerobic bacteria and archaea dominate the bottom waters and sediments. Remineralizing organic matter from the upper water column (mixolimnion) supports the microbes and the anaerobic oxidation of methane (AOM). Methane-oxidizing microbes, specifically anaerobic methane-oxidizing archaea (ANME), are coupled to sulfate-reducing bacteria in a microbial consortium (Boetius et al., 2000; Hinrichs and Boetius, 2002; Werne et al., 2004). ANME mediate methane levels in freshwater and soda lakes, and they account for approximately 90 % of methane consumed through AOM in modern oceanic systems (Egger et al., 2018). Rates of methane consumption differ by environment and type of ANME, with global freshwater systems ranging from 1 to 1 × 10⁵ $\mathrm{nmol}{\mathrm{CH}}_{\mathrm{4}}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ consumed (Martinez-Cruz et al., 2018). Although soda lakes have been less studied, consumption rates as high as 1.6 × 10⁴ $\mathrm{nmol}{\mathrm{CH}}_{\mathrm{4}}{\mathrm{L}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$ have been observed in freshwater Lake Kivu (Roland et al., 2018). Tracking AOM over geologic time periods is important because methane release from tropical wetlands was concomitant with the end of glacial conditions in Europe and is poorly constrained (DeMenocal et al., 2000; Riddell-Young et al., 2023). Additionally, large methane releases might have been partly responsible for the Permian–Triassic mass-extinction event (Berner, 2002).

Over geologic time, it is possible to gauge periods of increased methane oxidation, as shown by Zhang et al. (2011) in oceanic systems using a ratio of archaeal glycerol dialkyl glycerol tetraether (GDGT) lipids (de Rosa et al., 1977; Langworthy, 1977). The ratio, as described by Zhang et al. (2011), is known as the methane index (MI), and it uses GDGTs produced predominantly by euryarchaeal ANME. The MI has been employed to discern methanotrophy using the assumption that benthic methanotrophic Euryarchaeota preferentially produce GDGT-1, -2, and -3 and that the GDGTs crenarchaeol (cren) and crenarchaeol^′ (cren^′) are thought to come from Thaumarchaeota and Crenarchaeota, which are part of the TACK superphylum, typically found in the upper water column (Sinninghe Damsté et al., 2002; Pitcher et al., 2010; Zhang et al., 2011). Currently, the newly suggested names in the Genome Taxonomy Database for Thaumarchaeota and Crenarchaeota are Nitrososphaerota and Thermoproteota, respectively (Oren and Garrity, 2021; Rinke et al., 2021), which are used in this paper.

Kim and Zhang (2023) have shown a qualitative and quantitative relationship between the MI and methanotrophy in deep time, namely from the late Oligocene to the early Miocene. Kim and Zhang (2023) showed that the MI is applicable to AOM, with other biomarkers co-occurring in high-MI intervals representative of not only the Group I consortium of anaerobic methanotrophs (ANME) that produce GDGTs but also of Group II and Group III consortia (ANME-2 and ANME-3, respectively). To date, no studies have directly applied the MI to sediments in African soda lakes, despite evidence for AOM in modern soda lakes. Combined with MI values, other methane-related indices are used here to interpret methanogenesis and methanotrophy related to AOM. Previous studies have used GDGT-0 and GDGT-2 ratioed to the GDGT crenarchaeol value, which was originally thought to only be produced by mesophilic Thermoproteota (Blaga et al., 2009; Weijers et al., 2011). However, the optimum temperature for crenarchaeol production is closer to 40–45 °C (Zhang et al., 2007). Blaga et al. (2009) found that methanogens predominantly produced GDGT-0, whereas Weijers et al. (2011) showed that methanotrophic archaea predominantly produce GDGT-2.

Lake Magadi (Kenya) is a sulfate-limited lake; therefore, methanogenesis and methanotrophy may co-occur without suppression of the higher energy yield of sulfate reduction (Nijaguna, 2006; Sorokin et al., 2007; Deocampo and Renaut, 2016; Lameck et al., 2023). Here, we document evidence of methane cycling in Lake Magadi using archaeal isoprenoid GDGT lipid biomarkers. Environmental influences on archaeal community composition included precipitation/evaporation fluctuations and variations in hydrothermal activity, with the latter often related to contemporary tectonics. This study leverages four methane-related indices to understand methane cycling in recent and ancient lacustrine sediments: (1) the MI, (2) the %GDGT-0$/$crenarchaeol, (3) the %GDGT-2$/$crenarchaeol, and (4) the ratio of isoprenoid GDGTs $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ (hereafter $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$). Two distinct communities were found using a combination of the MI and ratios of GDGT-0 and GDGT-2, normalized to crenarchaeol. Intervals of high methanotrophy, as evidenced by MI and %GDGT-2$/$crenarchaeol, were related to an equally high proportion of methanogens, whereas the methane indices were low in periods when crenarchaeol was dominant.

Figure 1 Map of the drilling location of MAG-14-2A (yellow triangle) in Lake Magadi for the Hominin Sites and Paleolakes Drilling Project (HSPDP) (map credit: Esri 2024). Publisher's remark: please note that the above figure contains disputed territories.

2 Materials and methods

2.1 Study locations and sampling

Modern Lake Magadi is a seasonally flooded, saline alkaline pan composed of bedded trona ([Na₃(CO₃)(HCO₃) ⋅ 2H₂O]) located in the southern Kenya Rift near the border with Tanzania (Baker, 1958; Eugster, 1980). Its elevation is approximately 600–605 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ (above sea level), and it has a maximum depth during the rainy season of a few decimeters to ∼ 1 m (Fig. 1; Renaut and Owen, 2023). Although Lake Magadi is situated near the Equator, it lies in a rain shadow. Consequently, it has a large present-day moisture deficit (2400 mm evaporation versus 500 mm precipitation annually; Damnati and Taieb, 1995).

The modern alkaline lake is fed by ephemeral streams and alkaline hot springs (up to 86 °C at adjacent Nasikie Engida), distributed along faults around the shoreline (Baker, 1958; Crane, 1981; Allen et al., 1989; Renaut and Owen, 2023). Former high-level shorelines are preserved as coarse clastic sediments and (locally) as stromatolites around the lake. These vary in age and record lakes of different depths during the Quaternary. Outcrops of sediments relevant to this study and situated near Lake Magadi are the Oloronga Beds and the Green Beds. The chert-bearing Oloronga Beds in outcrop were lain down between ∼ 800 and 300 ka, with cores extending this back to 1 Ma (Owen et al., 2019). Green Beds outcrops include abundant chert and have been variably dated between 191 and 40 ka (Behr and Röhricht, 2000; Owen et al., 2019), with cores suggesting a range from 380 to 105 ka (Owen et al., 2019). The High Magadi Beds were deposited between ∼ 25 and 9 ka (Fairhead et al., 1972; Goetz and Hillaire-Marcel, 1992; Williamson et al., 1993; Behr and Röhricht, 2000; Owen et al., 2019; Reinhardt et al., 2019). Calcrete commonly caps the Oloronga Beds (Eugster, 1980), but fluvial erosion locally scoured those sediments before the Green Beds were laid down (Renaut and Owen, 2023).

Lake Magadi was cored as part of the Hominin Sites and Paleolakes Drilling Project (HSPDP) in June 2014 to further our understanding of the paleoenvironments in the East African Rift Valley, to better contextualize hominin remains and artifacts, and to understand possible environmental influences on hominin evolution and migration (Cohen et al., 2016). A 197.4 m core (MAG-14-2A) was drilled in the northern end of Lake Magadi (1°51^′5.76^′ S, 36°16^′45.84^′′ E; Owen et al., 2019). In total, 107.7 m of sediment was recovered, with an overall core recovery of 55.4 % (Cohen et al., 2016). Here, we use the age model from Owen et al. (2019). The core ranges from the modern trona surface (0 ka) to the Magadi Trachyte basement, dated to ∼ 1 Ma at the core site (Owen et al., 2019). Cores were sampled in 2016 during the initial core description at the Continental Scientific Drilling Facility (CSD, formerly LacCore) at the University of Minnesota, Minneapolis. Altogether, 61 samples, covering the period from 456 to 14.9 ka (Table S1 in the Supplement), were collected and freeze-dried from dark brown to black silty-clay intervals. Based on their color, these samples were expected to have a high total organic carbon that would yield the highest quantity of biomarkers for our study. The nonideal sampling strategy is due to poor core recovery.

Over the past million years, Lake Magadi has evolved through various phases: from swampy freshwater waterbodies to a large fresh to mildly saline lake that was continuously fed by rivers and groundwater to smaller hypersaline lakes bounded by the Magadi grabens that dried to trona pans partly fed by hot springs (Owen et al., 2019; Renaut and Owen, 2023). From 545 to 380 ka, the Magadi catchment progressively changed to a more arid condition, with the paleolake marked by abundant calcareous, organic-rich sediments (Owen et al., 2019). Periodic freshwater inundation occurred from 380 to 105 ka into a highly saline, alkaline lake that accumulated minor calcite and magnesium-rich calcite at the lake margins. Ash that fell into this waterbody reacted to form a variety of zeolites with anoxic, sulfate-rich bottom-water brines subjected to microbial sulfate reduction (Owen et al., 2019; Deocampo et al., 2022). The most recent phases of the lake (105 to 0 ka) were more evaporatively enriched, with abundant trona and minor nahcolite. Well-preserved diatoms in sediments deposited after ∼ 500 ka suggest very high aqueous silica in the paleolake in order to explain the preservation of their frustules under highly alkaline conditions, which may reflect strong evaporative concentration of silica-rich hydrothermal inflows (Owen et al., 2019).

2.2 GDGT and bulk organic preparation and analysis

2.2.1 Lipid extraction

To obtain a total lipid extract (TLE), 61 samples from Lake Magadi were freeze-dried and homogenized, and ca. 5–10 g of sediment was ultrasonically extracted with a 2:1 ratio of DCM:MeOH (where DCM represents dichloromethane). The TLE for each sample was treated with activated copper shot to remove elemental sulfur. The TLEs were then separated into three fractions (apolar, AP; polar one, P1; and polar two, P2) using activated alumina via short-column chromatography. The AP fraction was eluted with 4 mL of a 9:1 ratio of hexane (Hex):DCM ($v/v$), the P1 fraction was eluted with 4 mL of a 1:1 ratio of DCM:MeOH, and the P2 fraction was eluted with 4 mL of MeOH. The P1 fractions were dried down and redissolved in a 99:1 ratio of Hex:isopropanol (IPA) ($v/v$) and filtered through a 0.45 µm PTFE filter of 4 mm diameter prior to GDGT analysis.

2.2.2 Bulk organic δ¹³C_OM analysis

Samples were subsampled from the same intervals as organic biomarkers for bulk organic carbon isotope analysis. Powdered sediment samples were weighed in silver capsules, and carbonates were removed by adding 5 % HCl in 4 h increments. Samples were analyzed on a Costech elemental analyzer coupled to a Thermo Finnigan Delta V Plus isotope ratio monitoring mass spectrometer (IRMS). Samples are reported as per mil (‰) deviations from the Vienna Pee Dee Belemnite (VPDB) standard in conventional delta notation.

2.3 GDGT preparation and analysis

2.3.1 GDGT analysis

Polar samples from Lake Magadi were analyzed for core lipid isoprenoid glycerol dialkyl glycerol tetraethers (iso-GDGTs) at the University of Massachusetts Amherst on an Agilent 1260 series high-performance liquid chromatograph (HPLC; Fig. S1 in the Supplement) in tandem with an Agilent 6120 series single quadrupole mass selective detector (MSD). Compounds were ionized using atmospheric-pressure chemical ionization (APCI). The columns used for GDGT separation were a pre-column guard followed by two ultrahigh-performance liquid chromatography (UHPLC) silica columns (BEH HILIC, 2.1 × 150 mm, 1.7 µm, Waters) connected in series and kept at 30 °C. Elution solvents followed Hopmans et al. (2016) using a flow rate of 0.2 mL min⁻¹. Two solvent mixtures, hexane (A) and a 9:1 ratio of Hex:IPA (B), were eluted isocratically for 25 min with 18 % B, a linear gradient to 35 % B in 25 min, and a second linear gradient to 100 % B in 30 min.

2.3.2 GDGT indices

Several different ratios based on the relative abundance of different isoprenoid GDGTs have been developed to determine their source(s). The methane index (MI) is defined by Zhang et al. (2011) and is calculated as follows:

$$\begin{array}{}\text{(1)}& \text{MI}={\displaystyle \frac{\text{GDGT-1}+\text{GDGT-2}+\text{GDGT-3}}{\text{GDGT-1}+\text{GDGT-2}+\text{GDGT-3}+\text{Cren}+{\mathrm{Cren}}^{\prime }}}.\end{array}$$

MI values range between 0 and 1, with values > 0.5 considered to be derived from methanotrophic communities and values < 0.3 considered to be normal sedimentary conditions (Zhang et al., 2011). These proposed ranges from Zhang et al. (2011) were derived from GDGTs found in marine sediments, so the cutoff values for methanotrophy may differ in lacustrine sediments, particularly those in saline, alkaline environments.

The ratio of GDGT-2$/$crenarchaeol (%GDGT-2$/$cren) also indicates methanotrophy (values > 0.2), specifically methanotrophy associated with sulfate–methane transition zones (Weijers et al., 2011). These values were normalized and converted to percentages so that the numbers produced could be contextualized with the other indices used (Eq. 2). As a result, %GDGT-2$/$cren contributions greater than 33 % will be considered methanotrophic signals.

$$\begin{array}{}\text{(2)}& \mathit{\%}\text{GDGT-2}/\text{cren}={\displaystyle \frac{\left[\text{GDGT-2}\right]}{\left[\text{GDGT-2}\right]+\left[\text{Cren}\right]}}\times \mathrm{100}\end{array}$$

Methanogenic inputs are calculated similarly to Eq. (2) above using GDGT-0 in place of GDGT-2. Blaga et al. (2009) found that values of GDGT-0$/$(GDGT-0 + cren) > 2 are associated with methanogenic archaeal communities in a study of freshwater lakes. Similarly, in a study of Eocene marine sediments, Inglis et al. (2015) normalized the equation and converted it to a percentage, a convention that we follow (Eq. 3). They found that contributions from methanogens were indicated by values greater than 67 %.

$$\begin{array}{}\text{(3)}& \mathit{\%}\text{GDGT-0}/\text{cren}={\displaystyle \frac{\left[\text{GDGT-0}\right]}{\left[\text{GDGT-0}\right]+\left[\text{Cren}\right]}}\times \mathrm{100}\end{array}$$

The GDGT-2$/$GDGT-3 ($\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) index is used here to describe both mesophilic environments and environments with high MI values. Rattanasriampaipong et al. (2022) found that differences in $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ are linked to distinct archaeal communities, whereby low values of $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ (ca. 0.55) are observed in thermophilic cultures, whereas elevated values are indicative of hot-spring mats (ca. 1.00), shallow aerobic ammonia-oxidizing archaea (AOA; ca. 1.16), or archaea in suspended particulate matter (ca. 2.52). This is the same version described in Rattanasriampaipong et al. (2022) (Eq. 4).

$$\begin{array}{}\text{(4)}& \left[\mathrm{2}\right]/\left[\mathrm{3}\right]=\left[\text{GDGT-2}\right]/\left[\text{GDGT-3}\right]\end{array}$$

2.4 Bulk geochemistry

Bulk geochemical data and core descriptions from Owen et al. (2019, 2024b) were used to interpret hot-spring influences in the intervals of focus in the sediment core (designated as intervals 1, 3, and 5 in Sects. 3 and 4). They attributed rare Earth element (REE) anomalies to increased lake alkalinity, which reflected increased evaporation and the development of highly saline, alkaline lakes and possibly increased hydrothermal/fluvial inflow ratios. All statistical analyses were performed using the GraphPad Prism 10© software (https://www.graphpad.com/, GraphPad Prism version 10.3.1, GraphPad, 2025). Only the necessary data to determine relationships between the bulk geochemistry of REEs (La, Ce, Nd, Sm, Eu, Yb, and Lu) and methane indices (MI, %GDGT-0$/$cren, %GDGT-2$/$cren, and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) were imported. The data were tested for normality via the built-in “Normality (and lognormality) tests” function in the GraphPad software. Tests yielded lognormal distributions of each dataset and found the data to be non-normally distributed. Additionally, pyrite cubes were visually assessed without the use of a microscope during the initial core description and during organic extractions; however, no smear slides were collected or assessed.

2.4.1 Principal component analysis (PCA)

For the PCA, the imported data were analyzed using the built-in PCA function in the GraphPad software. The data were standardized, which scaled the data to a mean of 0 and a standard deviation of 1. The principal components were selected based on their eigenvalues using the Kaiser rule, which selects eigenvalues greater than 1.0. Principal components 1 and 2 explained 58.0 % of the variance in the data.

2.4.2 Correlation matrix

The correlation matrix was performed using the built-in function in the GraphPad software. Because the data were non-normally distributed, the nonparametric Spearman correlation was chosen over the Pearson correlation. An r value was computed for every pair of Y datasets using the default two-tailed option at a 95 % confidence interval.

3 Results

3.1 GDGT lipid variability

Samples are split into six intervals (1–6) based partly on their fractional abundances of GDGT-0 and cren as well as their MI values: (1) 35.67 to 32.61 m, ca. 17.7 to 14.9 ka; (2) 67.81 to 43.51 m, ca. 129 to 38.9 ka; (3) 86.06 to 70.78 m, ca. 197 to 149 ka; (4) 96.38 to 94.91 m, ca. 318 to 315 ka; (5) 104.10 to 103.16 m, ca. 324 to 323 ka; and (6) 130.21 to 119.64 m, ca. 456 to 391 ka (Table S1).

Figure 2 Downcore variations in Lake Magadi of the (a) %GDGT-0$/$cren, (b) MI, (c) %GDGT-2$/$cren, (d) bulk δ¹³C_OM, and (e) the GDGT-2$/$GDGT-3 ($\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) values from ca. 14.9 to 456 ka. Sections 1, 3, and 5 are low-MI intervals, outlined in blue, and sections 2, 4, and 6 are high-MI intervals, outlined in yellow. Bands at the top of the graph indicate the inferred (via Renaut and Owen, 2023) lake levels and major inputs, with dark blue indicating a perennial lake, brown indicating a sand or mud flat, and purple indicating tephra. Dotted lines in each section denote the cutoff points for methane-related indices: MI (> 0.5), %GDGT-2$/$cren (> 33 %), and %GDGT-0$/$cren (> 67 %). See Sect. 2.3.2 for more details. Note the breaks in the x-axis scale.

In each of the intervals of the core, MI, %GDGT-0$/$cren, and %GDGT-2$/$cren values oscillate between high and low values, changing abruptly from one interval to the next (Table S1). The methanotrophic (%GDGT-2$/$cren) and the methanogenic (%GDGT-0$/$cren) indices track similarly to MI values; that is, when values of MI are high, so are the other two indices. It should be noted that there are some large gaps in sampling between intervals in the core due to our sampling regime (i.e., targeting intervals with high apparent organic matter based on darker silty matrix). Interval 1 is characterized by a higher proportion of cren and lower overall index values. The %GDGT-0$/$cren index averages 36.3 % (± 0.09 %) in this interval, while the %GDGT-2$/$cren index averages 10.6 % (± 0.05 %; Fig. 2). The MI in this interval is correspondingly low, with an average of 0.25 (± 0.07), well below the MI = 0.5 cutoff range for methanotroph-impacted communities. As such, this interval could be used for the $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index: values averaged 2.1 (± 1.02). Interval 2 has much higher values for each of these indices, where the average %GDGT-0$/$cren = 99.3 % (± 0 %) and the average %GDGT-2$/$cren = 93.6 % (± 0.04 %). MI values in Interval 2 are also high, with an average of 0.96 (± 0.02). Of note, there is a large gap where no measurements were taken from 43.55 to 46.68 m (∼ 9.7 kyr) as well as from 50.36 to 58.74 m (∼ 32.5 kyr). Interval 3 averages for %GDGT-0$/$cren and %GDGT-2$/$cren are 54.3 % (± 0.27 %) and 20.4 % (± 0.26 %), respectively. However, there is one anomalously high value at 77.32 m, with %GDGT-0$/$cren and %GDGT-2$/$cren values at 99.6 % and 93.8 %, respectively, and an MI = 0.96. Excluding the high index value, the averages were lowered to 48.6 % (± 0.22 %) and 11.2 % (± 0.05 %) for the respective %GDGT-0$/$cren and %GDGT-2$/$cren values, and the MI average was lowered from 0.33 (± 0.23) to 0.26 (± 0.06). With the exclusion of 77.32 m, the $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index averaged 1.5 ( ± 0.80) in this interval, which is lower than in Interval 1. Interval 4 is characterized by high index values, with a similarly abrupt shift from low values. Averages of the %GDGT-0$/$cren and %GDGT-2$/$cren are 98.1 % (± 0.04 %) and 88.2 % (± 7.41 %), respectively, and an average MI value of 0.92 (± 0.04) is observed; these average index values are similarly high to those in Interval 2. Interval 5 is a shift to lower overall index values, with averages of %GDGT-0$/$cren and %GDGT-2$/$cren at 40.1 % (± 0.17 %) and 9.2 % (± 0.03 %), respectively, and an average MI of 0.22 ( ± 0.05). Finally, Interval 6 shows a period in the core with high index values throughout. Averages of %GDGT-0$/$cren, %GDGT-2$/$cren, and MI are 97.6 % (± 0.03 %), 89.4 % (± 0.08 %), and 0.95 (± 0.05), respectively.

3.2 Bulk δ¹³C_OM values

Bulk δ¹³C_OM values follow a similar pattern to the indices described in Sect. 3.1; that is, the values oscillate between high and low values between intervals. Samples in Interval 1 ranged from −21.9 ‰ to −16.8 ‰ and had an average δ¹³C_OM value of −18.4 ‰ with respect to VPDB. Interval 2 samples had the most ¹³C-depleted values in all sampled intervals, ranging from −89.4 ‰ to −24.7 ‰ with an average of −35.1 ‰; excluding the three outlier values (−48.1 ‰, −64.2 ‰, and −89.4 ‰), the Interval 2 average was −28.2 ‰. In Interval 3, the δ¹³C_OM had a narrower range from −24.4 ‰ to −21.4 ‰ and an average of −22.5 ‰. A lighter signal from Interval 4 yielded a narrow range of values from −27.0 ‰ to −25.4 ‰ and an average of −26.0 ‰. Interval 5 had slightly heavier values ranging from −25.0 ‰ to −18.1 ‰ with an average of −22.1 ‰. Lastly, Interval 6 had depleted δ¹³C_OM values similar to intervals 2 and 4, with a range of −28.2 ‰ to −22.1 ‰ and an average of −25.2 ‰. Analytical reproducibility of duplicate runs was better than ± 0.15 ‰ VPDB.

Figure 3 Downcore plot for Lake Magadi of (a) $\mathrm{Ca}/\mathrm{Na}$, (b) %REE abundance, and (c) MI values from ca. 14.9 to 456 ka. Sections 1, 3, and 5 reflect low-MI intervals, outlined in blue, and sections 2, 4, and 6 reflect high-MI intervals, outlined in yellow. Bands at the top of the graph indicate the inferred (via Renaut and Owen, 2023) lake levels and major inputs, with dark blue indicating a perennial lake, brown indicating a sand or mud flat, and purple indicating tephra. The dotted line on the MI plot (c) denotes the cutoff point > 0.5 for values significantly affected by methane-cycling archaea. Note the breaks in the x-axis scale and in the y-axis scale for the δ¹³C_OM. REE values are from Owen et al. (2019).

3.3 Bulk geochemistry

Both a PCA and correlation matrix were performed using the MI, $\mathrm{Ca}/\mathrm{Na}$, %GDGT-0$/$cren, %GDGT-2$/$cren, and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ to compare to the REEs La, Ce, Nd, Sm, Eu, Tb, Yb, and Lu (Fig. 3). Increased values of REEs are characteristic of sodic systems influenced by hydrothermal springs, such as Mono Lake in California and this system (Johannesson and Lyons, 1994; Owen et al., 2019). Additionally, $\mathrm{Ca}/\mathrm{Na}$ is a proxy for the proportion of hydrothermal inflow to Lake Magadi such that a higher $\mathrm{Ca}/\mathrm{Na}$ indicates more freshwater input, while a lower $\mathrm{Ca}/\mathrm{Na}$ indicates a higher proportion of hydrothermal water.

Figure 4 The (a) PCA and (b) Spearman correlation matrix showing the relationship between methane-related indices (MI and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) and REEs (La, Ce, Nd, Sm, Eu, Tb, Yb, and Lu) in the sampled intervals of the core. A negative relationship is seen between the methane indices and REEs, as shown by opposing eigenvectors on the PCA (a) and negative r values on the correlation matrix (b). REE values are from Owen et al. (2019).

A PCA (Fig. 4a) and nonparametric Spearman correlation matrix (Fig. 4b) were performed to quantify the relationship between REEs, MI, $\mathrm{Ca}/\mathrm{Na}$, and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$. The PCA showed that $\mathrm{Ca}/\mathrm{Na}$ loaded positively on PC1 and PC2 and that each of the methane indices loaded positively on PC1 and negatively on PC2. The REE La loaded positively on PC2 and negatively on PC1, while the REEs Ce, Nd, Sm, Eu, Yb, and Lu loaded negatively on PC 1 and PC 2. This indicates a negative relationship between the negatively loaded REEs and a high $\mathrm{Ca}/\mathrm{Na}$. Similarly, the correlation matrix of REEs and methane indices showed a negative relationship between each index and REE, except for the relationship of $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ and Nd (r = 0.02) and %GDGT-2$/$crem (r = −0.04), which showed no observable linear trend. The REEs and methane indices did not load on the same PC axis, showing that there was also not a nonlinear trend associated with the REEs and methane indices.

4 Discussion

4.1 Lake Magadi archaeal community shifts

The abrupt changes in isoprenoid GDGT-based indices in the sediment record of Lake Magadi indicate shifts in the archaeal communities present (Fig. 2). Shifts between two distinct communities were inferred using a combination of the methane index (MI) and ratios of GDGT-0 and GDGT-2 normalized to crenarchaeol (Eqs. 2 and 3). We denote these shifts as either “MI-on periods”, characterized by MI > 0.5 during intervals 2, 4, and 6, or “MI-off periods”, characterized by MI < 0.5 during intervals 1, 3, and 5. Oscillations between these two states are discussed in detail in the following sections.

Figure 5 Ternary plot of crenarchaeol, GDGT-0, and GDGT-2, which are used to calculate the methane indices. Samples are split by both their interval (denoted by their shape) and whether they are from a high-MI (yellow) or low-MI (blue) interval. Higher proportions of GDGT-0 indicate methanogenic inputs, higher proportions of GDGT-2 indicate methanotrophy, and higher proportions of crenarchaeol indicate more mesophilic conditions influenced by hot springs.

4.1.1 MI-on periods

In Lake Magadi, during the MI-on periods (Fig. 2; intervals 2, 4, and 6), the MI is persistently greater than 0.83 and displays more ¹³C-depleted δ¹³C_OM values compared to MI-off periods, indicating periods of enhanced methane cycling. AOM is a likely mode of methane cycling in Interval 6 and in parts of Interval 2, as sulfate-reducing bacteria (SRB) and AOM archaea live in a consortium together at the sulfate–methane transition zone, or SMTZ (Boetius et al., 2000; Hinrichs and Boetius, 2002; Werne et al., 2004). Thus, in intervals of the Magadi core where a SMTZ is suspected, such as in parts of Interval 2 and most of intervals 4 and 6, there should be an increase in indices related to methanotrophy, such as a high MI and %GDGT-2$/$cren (Weijers et al., 2011). Additionally, whereas methanogens and methanotrophs appear to be present in a consortium based on both the methane indices and bulk δ¹³C_OM, the majority of the contributions are coming from methanogens, as seen in the ternary plot in Fig. 5. This may seem counterintuitive, as the MI has been typically used to describe samples exhibiting a high predilection towards methanotrophy, but a high MI value does not necessarily exclude methanogenesis and, conversely, neither does a low MI; rather, a low MI value suggests a predominance of Thermoproteota over Euryarchaeota (Zhang et al., 2011). High %GDGT-0$/$cren and %GDGT-2$/$cren index values in intervals 2, 4, and 6 (Fig. 2) show that methanogenesis is co-occurring with AOM. The $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index is also useful for understanding the proportion of methanotrophs in sediments, even in intervals with high MI values like those discussed herein (Table S1; Fig. 2). Values of the GDGT $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ ratio track nearly identically to the MI values (Fig. 2), indicating that the MI is influenced by GDGT-2, which is characteristic of methanotrophs (Pancost et al., 2001; Schouten et al., 2003; Zhang et al., 2011).

Typically, methanogenesis in sulfate-rich systems is suppressed in favor of sulfate reduction caused by competition for both H₂ and organic substrates (Fazi et al., 2021; Sorokin et al., 2015). However, reports of methanogenesis co-occurring with SRB have been noted when methanogens are using noncompetitive substrates, such as methanol, or when sulfate levels are low (Oremland et al., 1982; Giani et al., 1984; Hoehler et al., 2001; Bebout et al., 2004; Arp et al., 2008, 2012; Jahnke et al., 2008; Smith et al., 2008; Robertson et al., 2009). Furthermore, pyrite cubes are common and scattered through the intervals where high index values are observed, indicating that there was a substrate for SRB, although it may have been in low concentration (Table S1). Thus, the combined evidence of pyrite in intervals with high GDGT-based indices (e.g., MI, %GDGT-0$/$cren, and %GDGT-2$/$cren $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) indicates the presence of a SMTZ that supports AOM with the co-occurrence of methanogenesis.

Looking at Fig. 5, it appears that GDGT-0 is the dominant GDGT compared to GDGT-2 and cren, indicating that this interval is likely methanogen-dominant rather than ANME-dominant. Interval 2 (Figs. 2 and 5) of the Magadi core appears to be more influenced by methanogenesis than AOM, as seen by a high %GDGT-0$/$cren signal accompanied by a high %GDGT-2$/$cren signal, high $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ ratios, and a more ¹³C-depleted bulk δ¹³C_OM signal (average = −35.1 ‰; median = −28.6 ‰). Values of bulk δ¹³C_OM are similarly ¹³C-depleted in AOM-dominant euryarchaeotal systems such as active mud volcanoes (∼ −27 ‰; ANME-1), a Danish freshwater lake (average ∼ −29.7 ‰; ANME-2), and the Sea of Galilee (∼ −30 ‰; ANME-2; Lee et al., 2018; Norði et al., 2013; Sivan et al., 2011). At points where the bulk δ¹³C_OM values are at their lowest (e.g., −89 ‰), they are accompanied by a lower %GDGT-2$/$cren at ca. 95 % and an elevated %GDGT-0$/$cren at > 99.5 %. This aligns with the literature, as Summons et al. (1998) reported δ¹³C_OM values of between −53.4 ‰ and −48.7 ‰ in the total lipid extract of methylotrophic methanogens using noncompetitive substrates in anoxic hypersaline environments. Furthermore, as these waters are typically sulfate-limited, acetoclastic and/or hydrogenotrophic methanogenesis is likely dominant when evidence of SRB is lacking (i.e., pyrite). Zhuang et al. (2016) performed compound-specific isotope analysis on several archaeol compounds from the Orca Basin and found archaeol and hydroxyarchaeol using H₂ or CO₂ (diagnostic of methanogens and methanotrophs) were relatively depleted (ca. −80 ‰ to −60 ‰) compared to the bulk ¹³C_OM (ca. −22 ‰). Zhuang et al. (2016) concluded that acetoclastic and/or hydrogenotrophic methanogenesis was unlikely due to high ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ concentrations in the Orca Basin, which may be the case in Lake Magadi. In Interval 2, there is no evidence of visible pyrite cubes, and we did not have a priority at the time of sampling to check a thin section of each sample for smaller pyrite aggregates. This indicates that other euryarchaeotal communities may have different forms of AOM occurring in the sediments. These other forms of AOM include nitrate/nitrite reduction and iron coupled to AOM (in `t Zandt et al., 2018). This is further bolstered by the evidence outlined by Kim and Zhang (2023) that not only quantitatively linked AOM to high MI values but also to non-Group I ANME Euryarchaea, as other non-GDGT-producing ANME (e.g., ANME-2 and ANME-3) were shown to co-exist with Group I ANME. In the intervals that are missing pyrite (i.e., most of Interval 2; see Table S1; Ferland, 2017), the pyrite may have either been too small to see with the naked eye or the excess H₂S could have been incorporated into the kerogen by reacting with labile organic matter. From 59.40 to 58.80 m, values of the bulk δ¹³C_OM dip as low as −89.4 ‰ (Table S1; Fig. 2), which aligns well with methanogenic archaeal biomass (Norði et al., 2013). However, as discussed above, there is likely acetoclastic and/or hydrogenotrophic methanogenesis co-occurring in these high-index intervals; this is likely the dominant process where sulfate-dependent AOM is absent, and the sulfate-dependent AOM is likely replaced by coupling to either nitrate/nitrite or iron reduction.

Samples in Interval 4 (Table S1) of the Magadi core have high index values but no evidence of sulfate-dependent methanotrophy except for high MI values. Thus, this interval is interpreted as being methanogenic rather than methanotrophic. The abundance of pyrite in the four samples with low MI values (Table S1; 104.10 to 103.16 m) indicates sulfate reduction not linked to AOM. This is not observed in any other level of the core, and a hypothesized series of reactions, which may be linked to an abundance of SRB, anaerobic ammonium-oxidizing (anammox) bacteria, and Thermoproteota (AOA) in the overlying water column, is described below. Due to periodic influxes of fresh water into Magadi, in addition to a permanent meromixis present in almost all samples in this study post 380 ka, the water column would have been oxic in the upper portion and anoxic below the chemocline. Freshwater pulses would have also brought nutrients to the lake, such as ammonia (${\mathrm{NH}}_{\mathrm{4}}^{+}$) and sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$). The oxic portion of the water column would have supported microaerophilic AOA that oxidize ${\mathrm{NH}}_{\mathrm{4}}^{+}$ to nitrite (${\mathrm{NO}}_{\mathrm{2}}^{-}$), which is then transported to the anoxic part of the water column (Straka et al., 2019). Here, anammox bacteria use excess ${\mathrm{NH}}_{\mathrm{4}}^{+}$ and ${\mathrm{NO}}_{\mathrm{2}}^{-}$ from the AOA and convert these to N₂. Excess ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ is simultaneously being used by SRB, creating HS⁻ that is reacting with iron species in the sediments and being buried as pyrite. Ladderane lipids characteristic of marine annamox bacteria (Jetten et al., 2009) were not studied in Magadi sediments. However, there is both 16S rRNA and lipid evidence for the production of ladderanes in hot springs in the western United States, suggesting that annamox bacteria can persist in hot-spring environments (Jaeschke et al., 2009). Additionally, Kambura et al. (2016) found evidence for Planctomycetes in both microbial-mat and water samples surrounding the hot springs of Lake Magadi, lending credence to the hypothesis of AOA persisting in Lake Magadi. Without other lines of evidence, however, these are hypothetical reactions for explaining excess pyrite in the sediments without accompanying MI values. Nonetheless, this explanation has merit because of the high relative abundance of both crenarchaeol and cren^′.

In nearly all of Interval 6 (Table S1), there is evidence for a higher proportion of methanotrophic archaea from 128.74 to 119.64 m (increased %GDGT-2$/$cren and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$) and methanogenesis in the intervals from 130.21 to 129.77 m (higher %GDGT-0$/$cren compared to %GDGT-2$/$cren; Table S1). Samples from 123.43 to 119.64 (Table S1) are of note because the $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ values are lower than the MI values, whereas MI and $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ values align nearly 1:1 in every other instance. This is likely due to GDGT-2 not being the dominant control of the MI. While both %GDGT-0$/$cren and %GDGT-2$/$cren are equally high, there may be other factors in the water column exporting GDGT-2 to the sediments, possibly from deep-dwelling Group I.1b Thermoproteota, although this is unlikely due to limitations of depth (Taylor et al., 2013). The 656 m paleoshoreline reported by Casanova (1986, 1987) and Casanova and Hillaire-Marcel (1987) would imply a maximum water depth of ∼ 50 m during the Late Pleistocene (African humid period – AHP) based on present topography. However, sedimentary evidence for such a high paleoshoreline is not seen throughout the Magadi Basin. Earlier water depths are also unclear, as accommodation space was changing as the axial rift developed with faulting and subsidence (Owen et al., 2024a). This is not deep enough (> 1 km) to support Group I.1b Thermoproteota per the constraints outlined in Taylor et al. (2013).

4.1.2 MI-off periods

In intervals characterized by low MI, %GDGT-0$/$cren, and %GDGT-2$/$cren values (MI-off intervals are odd numbers highlighted in blue in Fig. 2), the δ¹³C_OM values are ¹³C-enriched relative to those intervals characterized by higher index values (MI-on; Fig. 2). As the methane-cycling indices (%GDGT-0$/$cren and %GDGT-2$/$cren) are both predominantly influenced by the availability of crenarchaeol, MI-off periods are marked by increased production of crenarchaeol. Typically, crenarchaeol is produced in open-ocean systems, freshwater lakes, and soils by the mesotrophic aerobic ammonium-oxidizing phylum Nitrososphaerota. However, it can also be found in other environments, such as in hot-spring mats made by Thermoproteota (Pearson et al., 2004, 2008; Schouten et al., 2013). As Thermoproteota prokaryotes require oxygen to oxidize ammonium to nitrate, the increased presence of crenarchaeol in the MI-off intervals suggests periods when conditions were more oxic, at least in the upper water column. The increase in crenarchaeol and the low $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index values suggest that more Thermoproteota prokaryotes were present in Lake Magadi in those periods. As mentioned in Sect. 4.1.1, three groups of AOA are of interest for interpreting which archaeal groups are found in low-index intervals of Lake Magadi. Averages of $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ from the global dataset in Rattanasriampaipong et al. (2022) are as follows: hot-spring mats (mean = 1.00), shallow AOA cultures (mean = 1.16), and shallow core tops (mean = 2.64). Placing these on a continuum, we can approximate the environment from $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ averages in Magadi, although it should be noted that the shallow AOA and shallow core-top values in Rattanasriampaipong et al. (2022) are based on marine core tops, while the hot-spring mats are based on terrestrial hot springs like those observed around Lake Magadi (i.e., pH > 6.5).

Interval 1 captures a transition from a more arid East Africa to a wetter period at the onset of the AHP. During wetter periods, more allochthonous material is transported to the lake, which includes vegetal remains that impact the overall bulk δ¹³C_OM values. This allochthonous vegetation enriches the overall bulk δ¹³C_OM values more significantly than other intervals in the Magadi core. Average values of bulk δ¹³C_OM are −17.7 ‰ in Interval 1 and correspond to the δ¹³C_OM values of aquatic sedges mixed with a terrestrial signal of grassy woodland (Sikes, 1994). Pollen records in Lake Magadi indicate that a mixture of C₄ grassy woodlands and C₄ aquatic sedges was predominant in the landscape that surrounded Lake Magadi at this time (Muiruri et al., 2021). In agreement with the pollen record, the δ¹³C_OM values likely record a mixture of C₄ grasses and C₄ sedges. Similar δ¹³C values were reported in C₂₁ to C₃₃ n-alkanes in equatorial regions of Cameroon, ranging from −18.2 ‰ to −17.6 ‰ and recording the signals from C₄ grasses and sedges (Schwab et al., 2015). This all suggests that the bulk δ¹³C_OM signal is dominated by terrestrial biomass, unlike other sections of the core, and there does not appear to be a significant influence from the benthic microbial community (i.e., methane cyclers or SRB).

Values of the $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index average 2.1 in Interval 1, with some values as high as 3.74 and 4.63 at 33.28 and 33.03 m, respectively (Table S1). The higher values are closer to what is captured from deep oceanic suspended particulate matter (SPM) and deep-ocean core tops below the pycnocline, although caution should be used when comparing lacustrine and oceanic sediments (Rattanasriampaipong et al., 2022). The increase in %GDGT-0$/$cren (50.6 % and 54.3 %; Table S1) and the slightly increased MI values (0.37 and 0.41; Table S1) imply that these samples were deposited in a deeper lacustrine environment. Evidence for a deeper paleolake at ∼ 40 ft (∼ 17–18 m) above the modern lakeshore (Baker, 1958) is also observed in the High Magadi Formation (ca. 17.7 to 10.8 ka), indicating that there was fresh water flowing into the lake during the period of deposition in Interval 1, likely creating a fresher water cap on the meromictic Lake Magadi (Barker et al., 1991; Behr, 2002; Owen et al., 2019). However, excluding the high $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index values in Interval 1, the average is 1.6, which is closer to the hot-spring mats and shallow AOA cultures (Rattanasriampaipong et al., 2022). Likely, the higher $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ index values represent periods of increased methanogenesis occurring in the sediments, with AOA input from the upper water column likely induced by proportional increases in the amount of hydrothermal inflow to the lake (Sect. 4.2). Figure 4 shows that $\mathrm{Ca}/\mathrm{Na}$ is anticorrelated with REEs in both the PCA and correlation matrix. As the proportion of $\mathrm{Ca}/\mathrm{Na}$ decreases when REEs increase, we can say that, statistically, when it is drier (and thus proportionately more hydrothermally influenced) $\mathrm{Ca}/\mathrm{Na}$ decreases, REE values increase, and the methane indices are suppressed. In the periods of lower $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ values, the community is interpreted as being dominated by AOA and thermophilic AOA cultures (i.e., Thermoproteota; Rattanasriampaipong et al., 2022) and is further supported by high %cren and %cren^′. Kumar et al. (2019) described similarly low $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ values in the water column of Lake Malawi that are akin to values observed in Lake Magadi in both intervals 1 and 3. They found that values of a lower normalized $\left[\mathrm{2}\right]/[\mathrm{2}+\mathrm{3}]$, ranging from 0.55 to 0.59, in Lake Malawi were associated with the shallower Thermoproteota (Thaumarchaeota) Group I.1b. This is compared to higher values of $\left[\mathrm{2}\right]/[\mathrm{2}+\mathrm{3}]$ in the deeper-dwelling Thermoproteota Group I.1a, which means that most samples in Interval 1 are likely sourced from Group I.1b (Kumar et al., 2019). Kumar et al. (2019) concluded that Group I.1b Thermoproteota prokaryotes were contributing to the lower $\left[\mathrm{2}\right]/[\mathrm{2}+\mathrm{3}]$ values, while the deeper-dwelling Group I.1a Thermoproteota prokaryotes were more prevalent in aphotic portions of the water column (Kumar et al., 2019). The normalized $\left[\mathrm{2}\right]/[\mathrm{2}+\mathrm{3}]$ used by Kumar et al. (2019), with values ranging from 0.55 to 0.65, approximates values of $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ in the 1.30 to 1.65 range, as also described in this paper. More recently, Baxter et al. (2021) found that Thermoproteota I.1b prokaryotes are more prevalent in the upper oxygenated portion of the water column within the photic zone, as evidenced by a higher relative abundance of crenarchaeol and lower relative abundance of GDGT-2. Thus, our interpretations of thaumarcheotal AOA in Lake Magadi sediments are consistent with data from Baxter et al. (2021) and Kumar et al. (2019). This interpretation is consistent with Interval 1 being a period of proportionately higher freshwater input and ${\mathrm{HCO}}_{\mathrm{3}}^{-}$-rich hydrothermal input and a deeper lake overall, which would explain the accompanying increase in crenarchaeol.

Interval 3 $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ averages are lower overall (mean = 5.4; Table S1), with only one outlying high value (ca. 77.32 m at a value of 36.7). Excluding this high index value, the $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ average drops to 1.5, which is closer to what is observed in shallow AOA cultures and hot-spring mats. With most samples being closer to unity (i.e., $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ = 1.0), it is likely that hot springs had a greater influence on the community composition in these intervals. Samples that are closer to unity (70.78, 70.86, and 71.08 to 75.93 m) also have relatively ¹³C-enriched δ¹³C_OM values (mean = −21.8 ‰) compared to samples with a higher $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$. Average isotope values in Interval 3 are between oceanic hydrothermal vents (mean = −19.0 ‰) and terrestrial alkaline hot-spring systems such as the Bison Pool hot spring in Yellowstone National Park (mean = −24.9 ‰) (Shah et al., 2008; Schubotz et al., 2013). As elevated amounts of GDGT-2 (i.e., relative abundance > 45 %) are associated with Euryarchaeota and values in intervals 1, 3, and 5 are much lower than 45 % (Table S1), these intervals are likely dominated by thermoproteotal AOA (Pancost et al., 2001; Turich et al., 2007; Taylor et al., 2013). The archaeal community composition in intervals 1, 3, and 5 is independent of these external factors and is related to hydrothermal flows. This further supports the hot springs driving the lake archaeal community composition, as there was less overall precipitation and the thermoproteotal communities were more abundant during intervals 1, 3, and 5. Lastly, Interval 5, which only has four samples, has similarly low values of $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$ (average = 1.4) to intervals 1 and 3, which is likely indicative of thermoproteotal AOA cultures.

4.2 The influence of hot-spring runoff ratios on the archaeal methane cycles

Hydrothermal fluids in the basin are rich in carbonate and bicarbonate as well as Na⁺ ions, inferred to be a result of (1) the weathering and alteration of trachytic (silica-rich) basement rock (Jones et al., 1977; Allen et al., 1989) and (2) mantle-derived CO₂ discharged mainly along faults (Lee et al., 2017; Muirhead et al.). Ca²⁺ and Mg²⁺ are also very low (Deocampo and Renaut, 2016). Renaut and Owen (2023) note that hydrothermal waters become important contributors to lake recharge during arid phases when fluvial inflow declines. For example, Nasikie Engida, a small hypersaline lake northwest of Lake Magadi, is partly maintained today by hydrothermal inputs, with trona, nahcolite, and zeolites accumulating during dry periods when there is little or no fluvial inflow (De Cort et al., 2019; Renaut et al., 2021; Renaut and Owen, 2023). Magadi Core MAG14-2A lacks evidence of complete desiccation (e.g., mudcracks, soils, and calcrete), and the lake appears to have retained surface waters through multiple drought episodes during the last million years (Owen et al., 2019). In contrast, separate lakes in the neighboring Koora Basin (∼ 10 km to the east) have dried out many times, leading to soil formation in that basin (Owen et al., 2024a).

Owen et al. (2019) also noted that during periods when highly saline, alkaline water dominated at Magadi, ash was zeolitized, REE patterns developed significant anomalies, and pyrite developed in anoxic/euxinic bottom waters of a meromictic lake. Moreover, they noted the excellent preservation of diatom opaline silica in highly alkaline lakes after about 540 ka, which suggests that very high levels of silica in lake water preserved their frustules from dissolution. Hydrothermal springs at Magadi presently contain high silica, but their waters require further evaporation to achieve concentrations that would preserve diatom silica under very high pH conditions. High silica concentrations brought about by strong evaporation in a lake maintained by spring inflows (meteoric or hydrothermal) may partly explain the abundant chert through the Magadi sedimentary sequence.

Proportional increases in hot-spring water during periods of increased aridity would have favored the development of highly saline, alkaline waters with significant impacts on archaeal communities. Samples in Interval 6 represent a drier period when the lake area and volume had shrunk, partly due to tectonic influences (Owen et al., 2024a), and lake floor anoxia was prevalent (Owen et al., 2019). Between ∼ 380 and 105 ka (Interval 5 through mid-Interval 2; Fig. 2), the paleolake was frequently meromictic with a freshwater mixolimnion that supported freshwater planktonic diatoms, while the saline monimolimnion and lake floor waters favored the alteration of ash to a variety of zeolites (e.g., erionite, phillipsite, clinoptilolite, and analcime) (Owen et al., 2019). Similarly, from ∼ 105 to 0 ka (mid-Interval 2 through Interval 1; Fig. 2), low $\mathrm{Ca}/\mathrm{Na}$, increased Br, and abundant zeolite formation indicate saline conditions. After about 80 ka, tectonic adjustments and increasing aridity led to desiccation in the Koora Basin, suggesting that spring inflows were important in maintaining a hypersaline lake in the Magadi Basin, as they do today during dry seasons (Owen et al., 2019). $\mathrm{Ca}/\mathrm{Na}$ ratios have steadily decreased over time, although not at a 1:1 rate, which aligns with methane index values and sudden increases in these indices. Figure 4 shows that the $\mathrm{Ca}/\mathrm{Na}$ and methane indices are statistically different from one another in both the correlation matrix and PCA plot. In Fig. 4b, $\mathrm{Ca}/\mathrm{Na}$ is loaded positively on PC1 and PC2, while the methane indices ($\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$, MI, %GDGT-0$/$cren, and %GDGT-2$/$cren) are loaded positively on PC1 and negatively on PC2. Furthermore, REE data also appear to reflect our MI-off and MI-on periods, as the REEs are anticorrelated in the correlation matrix (Fig. 4b) and loaded on different PC axes (Fig. 4a). Consequently, we interpret the changes in salinity and alkalinity in the Magadi paleolakes as reflecting the impact of climate on spring runoff ratios into the Magadi Basin, which in turn have exerted significant impact on the archaeal communities through the last million years. Samples in the low-MI intervals (ca. 32.61–35.67, 70.78–75.93, and 103.2–104.1 m) likely reflect proportionally increased spring runoff ratios at Magadi caused by increased evaporation and decreased precipitation in the surrounding landscape.

5 Conclusions

Sediments in Lake Magadi track the environmentally driven changes in archaeal communities over the past ∼ 456 ka. Using the MI to track the predominantly archaeal inputs at Lake Magadi, we have observed sudden and distinct shifts between mixed communities of euryarchaeotal methanogens and methanotrophs transitioning to mesophilic AOA Thermoproteota communities and back again. This shift is driven, in part, by moisture balances in the East African Rift, with wetter conditions periodically causing freshwater floods into a saline lake to form a meromictic waterbody at Magadi, and with more archaea derived from the upper water column rather than the sediments, as evidenced by low MI, low $\left[\mathrm{2}\right]/\left[\mathrm{3}\right]$, and relatively ¹³C-enriched bulk δ¹³C_OM. Methane indices were typically higher during periods of reduced hydrothermal activity, indicating more euryarchaeal communities, whereas Thermoproteota communities thrived during periods of higher hydrothermal activity. This is a clear relationship between low MI values, spring runoff ratios, lake salinity, alkalinity, and the development of mesophilic Thermoproteota. This study is one of the first to examine methane cycling in a soda lake over geologic time and provides valuable insights into how variable these systems can be. Soda lakes are important ecosystems for methane cyclers and should be studied more closely so that we can improve the understanding of global methane contributions in the past and constrain methane sources in the future.