Large herbivore δ¹⁸O as a proxy for aridity in the South African winter and year-round rainfall zone

Abstract

This study explores patterning in δ¹⁸O values of tooth enamel in contemporary African herbivores from mainly C₃-dominated ecosystems. Evapotranspiration causes plants to lose H₂¹⁶O to a greater extent than H₂¹⁸O, leaving leaves enriched in ¹⁸O. In eastern Africa, ES species (evaporation-sensitive species: those obtaining water from food) tend to have more positive δ¹⁸O_enamel values than EI species (evaporation-insensitive species: those heavily dependent on drinking water); the magnitude of the difference increases with increasing aridity. We find the same pattern applies in the winter and year-round rainfall region of southern Africa, allowing us to use δ¹⁸O_enamel in fossil animals to examine paleo-aridity. We apply this approach to infer aridity at Quaternary fossil assemblages from present-day winter and year-round rainfall zones, including Elandsfontein (ca. 1–0.6 Ma), Hoedjiespunt (ca. 300–130 ka), and Nelson Bay Cave (23.5–3 ka). This analysis suggests that (1) at various times during the Pleistocene, Elandsfontein and Hoedjiespunt environments were wetter than last glacial maximum (LGM) to Holocene environments at Nelson Bay Cave (year-round rainfall zone); and (2) considered alongside other evidence from the year-round rainfall zone, wetter conditions across the Pleistocene–Holocene transition at Nelson Bay Cave suggests that climate changes at near-coastal sites may be out of phase with the adjacent interior.

Introduction

Oxygen isotope ratios of faunal bioapatite depend on those of food and drinking water, as well as on aspects of animal physiology (Bryant and Froelich, 1995; Kohn et al., 1996; Sponheimer and Lee-Thorp, 2001; Levin et al., 2006; Lee-Thorp, 2008). In eastern Africa, differences in δ¹⁸O between mammalian species that obtain most of their water from food and those that regularly drink water provide an index of aridity that is useful as a paleoenvironmental indicator (Levin et al., 2006; Blumenthal et al., 2017). Here, we explore the patterning in δ¹⁸O_enamel in contemporary large herbivores in the winter and year-round rainfall regions of southern Africa. We then compare these patterns with previously published δ¹⁸O_enamel datasets from Quaternary fossil assemblages in the present-day winter (Elandsfontein and Hoedjiespunt) and year-round (Nelson Bay Cave) rainfall zones. These results contribute to our understanding of Quaternary climate history in southern Africa and provide context for the later stages of human evolution in the region.

Since mammalian body temperature is approximately constant, the oxygen isotopic composition (¹⁸O/¹⁶O) of body tissues (here, tooth enamel) is directly related to δ¹⁸O of body water (Bryant and Froelich, 1995; Bryant et al., 1996; Lee-Thorp and Sponheimer, 2005; Lee-Thorp, 2008). Inputs of oxygen into the body include drinking water, atmospheric O₂, and water in food. Outputs include urine, feces, panting (loss of water vapor), perspiration, and respired CO₂. All inputs and outputs may be influenced by climate and/or physiology, as well as the animal's diet and behavior (Ayliffe et al., 1994; Bocherens et al., 1996; Sponheimer and Lee-Thorp, 1999; Hoppe, 2006; Murphy et al., 2007).

The δ¹⁸O of precipitation (which accounts for most source water for plants) is affected by the value of the oceanic source from which it originates, the amount of fractionation during evaporation, the fractionation that occurs while precipitation is developing, and finally, the trajectory of the associated air mass (McGuire and McDonnell, 2007). In plants, water is not isotopically fractionated when taken up and transported through the roots and suberized stem (Barbour et al., 2007). During evapotranspiration, preferential loss of H₂¹⁶O results in leaves becoming enriched in ¹⁸O (Marshall et al., 2007); this effect is greater at higher temperatures and/or lower precipitation (Barbour and Farquhar, 2000; Helliker and Ehleringer, 2000, 2002). Thus, herbivores that obtain a significant proportion of their water intake from leaves are expected to have higher δ¹⁸O values, especially in arid environments. In contrast, the δ¹⁸O values of herbivores that regularly drink surface water tend not to increase relative to the δ¹⁸O of meteoric water as aridity increases (Levin et al., 2006; Blumenthal et al., 2017).

Faith (2018) provided an overview of general patterns of feeding and drinking behavior among herbivores. Browsers and frugivores generally obtain a significant proportion (often all) of their moisture requirements from their food, lessening their dependence on surface water. Browsers and frugivores are generally enriched in ¹⁸O relative to grazers because leaf water tends to be ¹⁸O enriched relative to surface water (Kohn et al., 1996; Sponheimer and Lee-Thorp, 1999; Levin et al., 2006; Blumenthal et al., 2017). Grazers generally need to drink more frequently, especially if they consume dry grass. These patterns apply particularly to ruminants; non-ruminants (equids, rhinos, elephants, hippos) are more water dependent regardless of their diets (Faith, 2018). Differences can also be seen within the grazing subgroup. For example, grazers with behavioral and physiological adaptations to arid environments, most notably gemsbok (Oryx gazella) and oryx (Oryx beisa), tend to have limited dependence on surface water.

Drawing primarily from observations in eastern Africa, Levin et al. (2006) and Blumenthal et al. (2017) have demonstrated that ungulate δ¹⁸O_enamel can provide a proxy for aridity (strictly, water deficit, calculated as potential evapotranspiration minus mean annual precipitation at each site). Species in which the offset between δ¹⁸O_enamel and δ¹⁸O_{meteoric water} increased systematically with increasing water deficit were termed ‘evaporation sensitive’ or ES (e.g., Giraffidae, Hippotragini, and Tragelaphini), while those showed no increasing offset were termed ‘evaporation insensitive’ or EI (e.g., elephant, rhinos, suids, and equids) (Levin et al., 2006). The former are primarily ruminants that are independent of surface water (typically browsers and mixed feeders), and the latter include non-ruminants that are highly water dependent. The more arid the environment, the bigger the difference in δ¹⁸O_enamel between the two groups of animals.

To date, this approach has been applied only in low-latitude, summer-rainfall environments where the majority of grasses are C₄ (Levin et al., 2006; Blumenthal et al., 2017; Patterson et al., 2019). This study investigates whether a similar pattern applies in the contemporary winter and year-round rainfall region of southern Africa where the vast majority of plants, including grasses, are C₃. If so, it may offer a means of reconstructing aspects of paleoclimate at the many important archaeological and fossil sites in the region.

Methods

To provide a comparative baseline, we measured δ¹⁸O_carbonate on tooth enamel of 308 contemporary large herbivores from varied environments in the winter rainfall zone (WRZ) and year-round rainfall zone (YRZ) as well as in the immediately adjacent areas extending into the summer rainfall zone (SRZ) (sampling locations in Figure 1). These were included to allow for past shifts in the boundary of the WRZ. Herbivores sampled derive from the families Elephantidae, Rhinocerotidae, Equidae, Suidae, Hippopotamidae, Giraffidae, and Bovidae. They were obtained from areas of relatively undisturbed vegetation in the western and southern portions of South Africa. The majority were obtained from South African National Parks (SanParks) and the Western Cape Nature Conservation Board (Cape Nature). These animals died natural deaths and their remains were collected between January 2013 and July 2015. A handful of specimens were obtained from farms that preserve areas of natural vegetation. Animals known to have been translocated were avoided. Some sampled specimens came from the collection in the Department of Archaeology at the University of Cape Town. The complete list of samples and results is provided in Supplementary Tables 1 and 2.

Figure 1. Map of southern Africa indicating sample collection sites and the boundaries of the winter (W), summer (S) and year-round (YR) rainfall zones.

δ¹⁸O values for fossils from Quaternary fossil assemblages from Elandsfontein, Hoedjiespunt, and Nelson Bay Cave have been published previously (Luyt et al., 2000; Hare and Sealy, 2013; Lehmann et al., 2016; Sealy et al., 2020; details in those references). All data related to these sites have been combined and are available in Supplementary Table 3.

Enamel powder (~10 mg) was drilled from the entire height of each tooth, from the occlusal surface to the cervix. The sample therefore integrates the whole period of crown formation (except for enamel lost through wear), averaging seasonal variations. While attempts were made to sample molars, particularly M3, this was not always possible because we sampled what was available in the faunal collection due to curatorial constraints. However, the sampling strategy was to avoid heavily worn teeth. Enamel powder was pre-treated (as described by Lee-Thorp et al., 1997, and Sponheimer and Lee-Thorp, 1999), with modifications as in Luyt and Sealy (2018). Approximately 2 mg of pre-treated enamel was placed in a Thermo Finnigan Model II gas bench at 72°C and reacted with 100% H₃PO₄ for a minimum of 2.5 hours to convert carbonate to CO₂. Isotope ratios were measured on a Delta Plus XP isotope ratio mass spectrometer. Four aliquots of each of three standards NBS18 (δ¹⁸O_PDB = −23.20‰), NBS19 (δ¹⁸O_PDB = −2.20‰), and an internal standard Cavendish Marble (δ¹⁸O_PDB = −8.95‰) were analyzed along with 40 unknown samples in each run. For each unknown or standard, the mass spectrometer measured the evolved gas nine times, averaged to obtain the final value. A regression of the observed on the expected values for the three standards was used to correct values obtained for samples to the Vienna Peedee Belemnite (VPDB) scale. The standard deviation (square root of the variance) of repeated measurements of the standards is ≤ 0.2‰ (results of repeated measurements of standards are given in Supplementary Table 4). ¹⁸O/¹⁶O ratios are reported in the δ notation relative to the VPDB standard in parts per mil (‰).

Each sample was attributed to the biome in which it was found. Biomes are classified in terms of the principal vegetation type, which depends mostly on climate (Mucina and Rutherford, 2006). Different biomes sample different climatic conditions (see Table 1 for details of biome attributes). Environmental data were supplied by the South African Weather Service or taken from the literature. The variables used were mean annual precipitation (MAP), mean annual temperature (MAT), mean annual soil moisture stress (MASMS), mean annual potential evapotranspiration (MAPE), relative humidity (RH), summer aridity index (SAI), winter concentration of rainfall (WCR), and water deficit (WD) (see Table 2 for definitions and references).

Table 1. Meteorological variables and dominant vegetation types of biomes in southern Africa (Mucina and Rutherford, 2006). Refer to Table 2 for definitions of meteorological variables.

Table 2. Definitions of meteorological variables.

Water deficit, which is derived from MAP and MAPE, provides an indirect measure of annual water availability to plants (Gallego-Sala et al., 2010; Harrison et al., 2010); lower WD values indicate moister environments. One limitation of the study was the variable proximity between South African Weather Service stations and collection points (ranging from 5 to 90 km). Reliance on environmental variables at the scale of the biome/bioregion underestimates variation across smaller geographical scales, but since our purpose here is to reconstruct paleoenvironments, we consider this broad-brush approach reasonable.

Univariate regressions (ranked ANCOVA, for non-parametric data) were used to explore the relationship between δ¹⁸O and each meteorological variable. Spearman's rho correlations were used to measure the strength of association between each pair of variables. Correlation coefficients were considered significant at p < 0.05. Animals were grouped by biome, feeding preference (grazers, browsers, mixed feeders; after Skinner and Chimimba, 2005), and sensitivity to evaporation (ES or EI) (Table 4). The independent categorical variables are ES or EI. The dependent continuous variable is the δ¹⁸O value and the continuous covariates are the meteorological variables.

Classification of species as either ES or EI is central to development of the δ¹⁸O aridity index. To provide such assignments, we relied heavily on the digestive physiology and drinking requirements of different species (data from Hempson et al., 2015), which are key determinants of their sensitivity to evaporation (Faith, 2018). All ruminant ungulates that are highly water dependent, which in our dataset included only grazers (e.g., Connochaetes taurinus, Syncerus caffer, Redunca arundinum), are considered to be EI (Table 4). We also included all hindgut fermenters (non-ruminants) in the EI category, because they tend to be water dependent regardless of their diets (Faith, 2018) and they are consistently insensitive to evaporation (Blumenthal et al., 2017). Ruminant ungulates that are independent of surface water or that have low requirements for surface water, including various browsers (e.g., Sylvicapra grimmia, Tragelaphus scriptus, Tragelaphus oryx), a mixed feeder (Antidorcas marsupialis), and an arid-adapted grazer (Oryx gazella), were assigned to the ES category. These assignments are broadly supported by empirical observations from tropical African environments (Blumenthal et al., 2017).

It is important to note that obtaining local δ¹⁸O_water values was outside the scope of this study and therefore we did not calculate ‘enrichment factors’ (i.e., enrichment in ¹⁸O_enamel compared with local meteoric water) (Levin et al., 2006; Blumenthal et al., 2017). The climatic and environmental gradients in our study area are much smaller than those encompassed by Levin et al. (2006) and Blumenthal et al. (2017), which ranged from arid settings (e.g., Lake Turkana) to high-altitude environments where annual precipitation considerably exceeds potential evapotranspiration (e.g., the Aberdare Range). We assumed that δ¹⁸O_water remains approximately constant across each biome and based our interpretations on differences between ¹⁸O_enamel of ES and EI animals per biome. In support of this approach, we note that the weighted average for δ¹⁸O_rainfall in Cape Town for the period 1996–2008 was −3.3 ± 1.8‰ (Harris et al., 2010), while at Mossel Bay, approximately 330 km east of Cape Town, the figure for 2009–2012 was −2.7 ± 2.44‰ (Braun et al., 2017). At Elandsfontein, approximately 100 km north of Cape Town (Fig. 1), the mean δ¹⁸O value for spring, tap, and standing water is −1.7 ± 2.2‰ (n = 4) (Lehmann et al., 2016). All three localities are in the Fynbos Biome. The δ¹⁸O_rainfall values are very similar, even though Cape Town and Elandsfontein receive mainly winter rainfall while Mossel Bay receives year-round rainfall, with a greater proportion originating over the Indian rather than the Atlantic Ocean (Bradshaw and Cowling, 2014).

Note that these δ¹⁸O values for water are relative to the standard mean ocean water (SMOW), not the VPDB standard. They can therefore be compared with one another but need adjustment before comparison with values for fauna reported below. Formula for conversion of SMOW to VPDB is (also given in the Supplementary Table 3):

δ¹⁸O_vPDB = (δ¹⁸O_vSMOW − 30.91)/1.03091 (Sharp, 2006)

The fossil samples from Hoedjiespunt 1 that we report here come from the paleontological assemblages accumulated by brown hyaenas (Parahyaena brunnea) and excavated between 1993–1998 (Berger and Parkington, 1995; Stynder, 1997; Parkington, 2003; Parkington et al., 2004). Some samples are from the surface-collected “old Hoedjiespunt” (Klein, 1983) collections, which show species composition and surface modifications indicating that they too probably derive from the hyaena den (Cruz-Uribe, 1991). An infrared-stimulated luminescence (IRSL) date of 345 ± 31 ka was obtained from sediment surrounding the bones and teeth, and a U-series date of circa 300 ka was obtained for the calcrete capping the deposit (Berger and Parkington, 1995; Stynder, 1997; Stynder et al., 2001). No samples used in this paper derive from the overlying archaeological shell middens, recently dated by optically stimulated luminescence (OSL) to the last interglacial (Tribolo et al., 2022). The hyenas burrowed into sediments already in place, which means that the fossil assemblage is younger than the surrounding matrix and likely dates to between ca. 300 ka and 130 ka.

The samples reported in this study for Elandsfontein derive from the Elandsfontein Main (EFTM) surface collections (Klein et al., 2007) as well as more recently excavated collections from the West Coast Research Project (WCRP) (Braun et al., 2013). Comparisons of the fauna found at the site with well-dated eastern African faunas suggest an age of between 1.0 Ma and 600 ka (Klein et al., 2007), with more recent paleomagnetic evidence indicating a minimum age of 780 ka (Braun et al., 2013). The potential for considerable time-averaging at Elandsfontein poses potential challenges for paleoenvironmental reconstruction because the fossil assemblage may sample different periods of time and varied environmental conditions.

Nelson Bay Cave is situated in the Southern Cape, currently in the year-round rainfall zone. Excavations in the 1960s and 1970s yielded rich faunal assemblages dating to the last 23 kyr (Inskeep, 1987; Klein, 1972). δ¹⁸O values previously reported for this fauna (Sealy et al., 2020) have been aggregated into three time periods: last glacial maximum (LGM), the last glacial–interglacial transition (LGIT), and the Holocene (chronology from Loftus et al., 2016). During the LGM, the site was approximately 30 km from the coastline (Carr et al., 2019), while in the two latter time periods it was on the coast.

Relevant fossil δ¹⁸O datasets from the region are also available from Boomplaas Cave (Sealy et al., 2016) and Elands Bay Cave (Stowe and Sealy, 2016), but these are not considered here. Small samples of ES taxa at Boomplaas Cave preclude application of the δ¹⁸O aridity index. Additionally, at Elands Bay Cave, we have noticed discrepancies between the originally published faunal lists (Klein and Cruz-Uribe, 1983, 2016) and the reported taxa for which δ¹⁸O values have been measured (Stowe and Sealy, 2016). Pending resolution of these discrepancies, we are hesitant to examine δ¹⁸O values as an indicator of aridity at the site.

Results

Biome-level variations by feeding type

To investigate possible differences in δ¹⁸O by biome, herbivores were divided into groups based on diet preference. For browsers (Table 3, Fig. 2), δ¹⁸O varied significantly by biome (Kruskal–Wallis H (5) = 60.76, p < 0.001). The most-positive values occurred in the Savanna biome (median = 6.1‰, n = 21) and the most-negative values occurred in the Forest biome (median = −4.6‰, n = 21). Pairwise Wilcoxon tests demonstrate significant differences between the Savanna biome and the Albany Thicket (p < 0.001), Fynbos (p = 0.001), and Succulent Karoo (p = 0.002). Values for the Forest biome were significantly different (p < 0.001) from all other biomes except Albany Thicket (p = 0.066). Fynbos and Succulent Karoo biomes (those with the largest proportions of winter rainfall) displayed similar medians of 1.0‰ and 1.8‰, respectively (p = 0.150).

Table 3. Summary of δ¹⁸O by feeder type and biome.

Figure 2. δ¹⁸O (‰) of browsers, grazers, mixed feeders, and omnivorous ungulates by Biome. Horizontal line in each box indicates median, X indicates mean. Sample sizes are in Table 3.

For the grazing species (Table 3, Fig. 2), δ¹⁸O also varied significantly by biome (Kruskal–Wallis H (4) = 63.44, p < 0.001). The most-positive values derived, once again, from the Savanna biome (median = 4.1‰, n = 36) while the lowest median values were from the Albany Thicket (median = −0.1‰, n = 52) and Fynbos biomes (median = −0.5‰, n = 12), although the sample size in the latter was small. There were no grazers from the Forest biome. Values for the Savanna biome differed significantly from those for the Fynbos biome and the Albany Thicket (p < 0.001), while values for these last two overlapped. There was also a significant difference between the Albany Thicket and the Succulent Karoo (p < 0.001).

The mixed feeders consist only of springbok (Antidorcas marsupialis), which in our dataset occurs in three biomes (Nama Karoo, Succulent Karoo, and Savanna), all of which are arid (Table 3, Fig. 2). A Kruskal–Wallis test revealed that the δ¹⁸O distributions of animals from these different biomes are indistinguishable (H (2) = 1.02, p = 0.601), as expected given the uniformly arid conditions represented in the sample.

Biome-level variations by sensitivity to evaporation

Figure 3 plots the δ¹⁸O values for each of the ES and EI groups of animals from each biome (see classifications in Table 4). The biggest difference between ES and EI animals was 5.5‰ in the Nama Karoo biome (Table 5, Fig. 3), followed by 2.6‰ in the Savanna biome. In both the Succulent Karoo and Forest biomes, samples of ES and EI species are small and not significantly different from each other. The highly diverse Fynbos region requires a higher resolution study to understand it thoroughly. The overall pattern is that the difference between ES and EI animals is small in cooler and/or moister biomes and increases with increasing aridity. Interpretation is limited by variations in sample sizes and species compositions.

Figure 3. Differences between mean δ¹⁸O of evaporation-sensitive (ES) and evaporation-insensitive (EI) animals for each biome. Horizontal line in each box indicates median, X indicates mean.

Table 4. Summary of ES/EI classifications for modern samples.

Table 5. Median δ¹⁸O_enamel per biome for the evaporation-insensitive (EI) and evaporation-sensitive (ES) animal groupings and the difference (ES − EI).

Relationships between δ¹⁸O and meteorological variables

Relationships between δ¹⁸O and meteorological variables were investigated by exploring correlations and then constructing regression models. Table 6 provides the correlation coefficients and significance of the relationships between the meteorological variables and δ¹⁸O values by feeding type and evaporation sensitivity. A summary of the outcomes is provided in Table 7.

Table 6. Spearman's correlation coefficients for δ¹⁸O and meteorological variables from all locations. (a) By feeder type; (b) by ES/EI category. **Correlation is significant at the 0.01 level; *Correlation is significant at the 0.05 level.

Table 7. Summary outcome of regression models for δ¹⁸O_enamel against each meteorological factor. (a) All groups; (b) ES group; (c) EI group. B is the unstandardized regression coefficient (the change in the isotope value with a one unit change in meteorological factor); t is the test statistic; abbreviations for meteorological variables explained in Methods.

δ¹⁸O values for grazers and browsers were significantly correlated with all meteorological variables. Grazer δ¹⁸O had highest correlation coefficients with MAP, MASMS, WD, MI, and MAPE, while browser δ¹⁸O was highly correlated with MAP, MAPE, WD, and MI. Mixed feeders showed few significant correlations, with only MAP and MI being significant, although our observations pertain to a single species (Antidorcas marsupialis) from arid environmental contexts, which limits our ability to detect meaningful correlations.

The δ¹⁸O values for the ES and EI taxa are significantly correlated with all of the meteorological variables examined here, with the strongest correlations observed for MAP, MAPE, and WD (Table 6, Fig. 4). For two of these variables (WD and MAPE), we observed significant differences in the slope between ES and EI taxa (WD: F = 4.415, p = 0.036; MAPE: F = 4.652, p = 0.032; MAP: F = 2.6317, p = 0.106; MI: F = 3.6231, p = 0.0579), such that δ¹⁸O increases at a faster rate in ES taxa as aridity increases (Fig. 4). The result is that the difference in δ¹⁸O between ES and EI taxa becomes larger in more arid environments.

Figure 4. Scatter plot of herbivore δ¹⁸O against the best meteorological variables by EI and ES. (a) Scatter plot of herbivore δ¹⁸O against WD by EI and ES; (b) scatter plot of herbivore δ¹⁸O against MAP by EI and ES; (c) scatter plot of herbivore δ¹⁸O against MAPE by EI and ES.

Discussion

Our analyses demonstrate that ES animals are enriched in ¹⁸O compared with EI ones, both within particular biomes (Fig. 3) and across our entire sample (Fig. 4). This pattern is in accordance with previous studies, which have found that animals that consume evaporatively enriched leaves have more positive δ¹⁸O values (Levin et al., 2006; Lee-Thorp, 2008; Braun et al., 2013; Lehmann et al., 2016; Robinson et al., 2016; Blumenthal et al., 2017). We did not see significant differences between ES and EI animals in the Forest biome (Fig. 4), which is characterized by the lowest potential evapotranspiration (MAPE) in our sample. It is likely that in this context, the evaporative enrichment of ¹⁸O in leaves relative to drinking water is negligible. The considerable overlap between ES and EI taxa in the Succulent Karoo biome, appears to relate in part to the massive variation in δ¹⁸O among some of the ES taxa. For example, the common duiker (Sylvicapra grimmia) from a single locality (Anysberg) have δ¹⁸O values ranging from −8.10‰ to 6.66‰.

Comparison of browsers across biomes shows that the most-positive values occur in the Savanna biome and the most-negative values occur in the Forest biome. Samples collected from the Forest biome came from close to the coast and represent the most mesic conditions in our sample, whereas those from the Savanna biome came from much farther inland in a more arid setting (Fig. 1). The patterns can be explained by both gradients in δ¹⁸O of meteoric water and local environmental effects within each biome. The more-positive values from the Savanna could be due in part to the increasingly positive δ¹⁸O in precipitation with distance from the source of rain (the continental effect) (Dansgaard, 1964; McGuire and McDonnell, 2007). This pattern was also noted at archaeological sites Nelson Bay Cave (on the coast) and Boomplaas Cave (inland) where the mean δ¹⁸O values for browsers are reported to be −2.8‰ and 3.1‰, respectively (Sealy et al., 2016, 2020). Values for grazers follow the same pattern, with means of −5.1‰ at Nelson Bay Cave and 1.5‰ at Boomplaas (Sealy et al., 2016, 2020).

δ¹⁸O values for grazers show less variation across biomes than browsers, but a similar pattern is evident with the highest values in the Savanna biome. There are no grazers from the Forest biome, so no comparison can be made there. The single species of mixed feeder (the springbok Antidorcas marsupialis) had relatively positive δ¹⁸O values, with only a small range of variation. Dietary studies (summarized in Skinner and Chimimba, 2005) show that springbok mostly browse on evaporatively enriched leaves in the arid environments in which they live, grazing only when there is fresh, new grass after rain. Therefore, it is not surprising that positive δ¹⁸O values, clustering with browsers, were noted in the Nama Karoo and Succulent Karoo biomes. This tendency towards browsing was also shown by their δ¹³C values (Luyt et al., 2019). Relatively positive δ¹⁸O values for Antidorcas marsupialis (clustering with browsers rather than grazers) have also been reported from Equus Cave by Sponheimer and Lee-Thorp (1999).

The best predictors for herbivore δ¹⁸O seem to be WD, MAPE, and MAP, with R² (adj) values of 0.54, 0.51, and 0.51, respectively. Ayliffe and Chivas (1990), Levin et al. (2006), and Murphy et al. (2007) all used δ¹⁸O of bioapatite as a proxy for environmental variables such as water availability and rainfall. This is consistent with the Levin et al. (2006) aridity index, which found larger differences between the δ¹⁸O values of ES and EI animals in more arid environments (higher WD values). The data from this study indicate that a variation of the aridity index does indeed apply to the winter and year-round rainfall zones of southern Africa. The difference in δ¹⁸O between animals sensitive to evaporation levels and those insensitive to evaporation levels increases in more arid environments (Figs. 3 and 4).

This study shows that the difference between the δ¹⁸O values of water-dependent/evaporation-insensitive animals increases with increased aridity. The results are sufficiently promising that in the future it would be worth measuring the meteoric water values for different localities across the winter and year-round rainfall zones. A further avenue for study would be serial sampling of these teeth. A limitation of this study is that while bulk sampling of teeth is a valuable method, it has inherent limitations. One notable constraint is that bulk samples may not be representative of a full year. In equids, for example, tooth crowns mineralize over periods of up to two years (Hoppe et al., 2004), but crowns in smaller mammals develop over shorter times. Tooth wear introduces another confounding variable because it erases portions of the chronological record stored in the crown. By serially sampling teeth one could be more confident of obtaining values over the entire time the tooth was mineralizing.

δ¹⁸O values as a measure of Quaternary aridity in the winter and year-round rainfall zones

Our examination of present-day herbivores demonstrates that the difference in the δ¹⁸O values of ES taxa relative to EI taxa ($\varepsilon _{ES-EI}$) increases in more arid environments. Thus, we expect lower $\varepsilon _{ES-EI}$ in wetter environments, and greater $\varepsilon _{ES-EI}$ in more arid environments. To illustrate how this δ¹⁸O aridity index might be applied to the fossil record, we re-assessed published δ¹⁸O values for two fossil faunal assemblages from the winter rainfall zone (Elandsfontein and Hoedjiespunt) and one from the year-round rainfall zone (Nelson Bay Cave). ES and EI assignments for fossil taxa were based on those for the nearest living relatives (see Supplementary Table 3 for assignments).

To estimate uncertainty in the $\varepsilon _{ES-EI}$ values across the fossil assemblages and to facilitate comparisons between them, we implemented a bootstrap approach using the R software (v.4.3.1; R Core Team, 2023). For each assemblage: (1) we generated a bootstrap sample (i.e., sample with replacement) of δ¹⁸O values for both ES and EI taxa; (2) for each bootstrap sample, we calculated the difference between the mean δ¹⁸O of ES taxa and the mean δ¹⁸O of EI taxa ($\varepsilon _{ES-EI}$); and (3) we repeated the process over 5000 iterations to generate a sampling distribution of $\varepsilon _{ES-EI}$ values from which a 95% confidence interval (CI) could be calculated.

Results of this analysis are reported in Table 8 and Figure 5. The degree of uncertainty surrounding our estimate of $\varepsilon _{ES-EI}$ is comparable across assemblages, ranging from 1.3‰ (Nelson Bay Cave: Holocene) to 2.3‰ (Nelson Bay Cave: LGIT). The two winter rainfall zone assemblages have overlapping 95% CI for $\varepsilon _{ES-EI}$, such that they are indistinguishable from each other. Importantly, both sites provide a wetter signal (lower $\varepsilon _{ES-EI}$) than at Nelson Bay Cave. The 95% CI for Hoedjiespunt is lower than that for all time bins at Nelson Bay Cave (LGM, LGIT, and Holocene). At Elandsfontein, the 95% CI overlaps with the Holocene and LGM assemblages from Nelson Bay Cave, but it is significantly lower than the LGIT assemblage. At Nelson Bay Cave, there are detectable temporal changes, with the LGIT providing a significantly more arid signal than the Holocene. The 95% CI for the LGM assemblage from Nelson Bay Cave overlaps with both the LGIT and Holocene, such that we could not distinguish it from either of the other two time bins using the δ¹⁸O aridity index.

Table 8. The difference in the δ¹⁸O values of evaporation-sensitive (ES) taxa relative to evaporation-insensitive (EI) taxa for Elandsfontein, Hoedjiespunt, and Nelson Bay Cave (LGM = last glacial maximum; LGIT = last glacial-interglacial transition). 95% Confidence intervals are based on the bootstrap procedure described in the text.

Figure 5. Smoothed histograms comparing bootstrapped δ¹⁸O aridity index values (mean δ¹⁸O of ES taxa minus mean δ¹⁸O of EI taxa: $\varepsilon _{ES-EI}$) for Nelson Bay Cave (NBC), Hoedjiespunt (HDP), and Elandsfontein (EFT). More-positive values indicate more-arid conditions. LGM = last glacial maximum; LGIT = last glacial–interglacial transition; HOLO = Holocene.

There are several broader paleoenvironmental implications of these results. First, we note that the present-day environment at the two winter rainfall zone sites is considerably more arid than at Nelson Bay Cave. Thus, the observation of relatively wet conditions compared to Nelson Bay Cave implies substantial paleoenvironmental change in the winter rainfall zone at various times during the Pleistocene. The well-watered environments implied by the δ¹⁸O aridity index is consistent with multiple lines of evidence from Hoedjiespunt and Elandsfontein. Both sites include species that are associated with standing water, such as southern reedbuck (Redunca arundinum) or hippopotamus (Hippopotamus amphibious) (Klein, 1983; Klein et al., 2007). Additionally, ecometric analysis of the distribution of dental traits across the fossil herbivore communities indicates that Elandsfontein and Hoedjiespunt were both wetter than present-day environments in the winter and year-round rainfall zones of southernmost Africa (Faith et al., 2020). At Elandsfontein, the presence of wetland phytoliths (Mann, 2017) and spring deposits (Braun et al., 2013) provide further evidence for standing water.

Although there is ample evidence that the paleoenvironments sampled by the faunas from Hoedjiespunt and Elandsfontein were wetter than today, including large bodies of standing water, the precise mechanisms that drove these changes are less clear. Glacial phases of the Pleistocene generally were associated with wetter conditions in the winter rainfall zone (e.g., Chase and Meadows, 2007; Chase et al., 2017), whereby expanded Antarctic sea ice translates to a northward shift of the southern hemisphere westerlies, bringing more frequent winter storms to southern Africa. This is reinforced by downscaled LGM climate models that suggest considerable increases in winter rainfall zone precipitation (Engelbrecht et al., 2019), although nearby speleothems suggest that climatic changes during middle Pleistocene glacial–interglacial cycles (perhaps overlapping in age with Hoedjiespunt) may have been less pronounced (Braun et al., 2023). A complicating factor is that the cooler conditions thought to be associated with greater rainfall should also be linked to lower sea levels, which would contribute to a lower water table—yet there is clear evidence for standing water at both fossil sites. At Hoedjiespunt, the rarity of marine fauna is thought to signal a reduction in sea level (Stynder, 1997). Thus, to drive an expansion of wetlands, the amount of additional moisture in the environment would have to offset the water table declines linked to marine regression. It may be that cooler temperature, and the associated reduction in evapotranspiration, played an important role in facilitating greater moisture availability and promoting wetland environments (e.g., Faith et al., 2020).

In addition, our observations provide insight into Late Quaternary climate change on the southern Cape coast at Nelson Bay Cave (year-round rainfall zone). The overlapping 95% CI for both the LGM and Holocene implies relatively muted differences in moisture availability under LGM conditions. This is inconsistent with observations from inland sites in the year-round rainfall zone, which suggest considerable increases in moisture availability from the LGM to the Holocene (e.g., Chase et al., 2017; Faith et al., 2019, 2024). However, there is emerging evidence that climate changes on the southern Cape coast may be out of phase with sites in the adjacent interior, perhaps due to the localized influence of the warm Agulhas current on near-coastal environments (Chase and Quick, 2018). LGM climate models further suggest spatial variation in the year-round rainfall zone, with decreased rainfall on the southern Cape coast coupled with increased rainfall in the adjacent interior (Engelbrecht et al., 2019). With this in mind, our observation of a relatively dry LGIT followed by a wetter Holocene—which is the opposite signal observed at inland sites in the year-round rainfall zone (Chase et al., 2017; Faith et al., 2019, 2024)—supports the hypothesis of an anti-phase relationship between inland and coastal sites within the year-round rainfall zone (Chase and Quick, 2018).

Conclusion

Large herbivore δ¹⁸O can be used as a proxy for aridity in the South African winter and year-round rainfall zone. The approach taken in this study differs from that of Levin et al. (2006) and Blumenthal et al. (2016), in that it does not require knowledge of δ¹⁸O_{precipitation}, either through direct measurement or the use of aquatic fossil fauna (not always present in fossil assemblages) as a proxy. Instead, we used the difference in δ¹⁸O between evaporation-sensitive and evaporation-insensitive fauna to examine aridity. This approach allows for wider application in paleo contexts where δ¹⁸O_{precipitation} is unknown. A more intensive sampling program would clarify the patterns seen in the current dataset and might enable more precise paleoenvironmental reconstruction.

Re-analysis of δ¹⁸O values reported from Quaternary fossil assemblages from the winter and year-round rainfall zones provide two key insights: (1) at various times during the Pleistocene, the winter rainfall zone faunas (Elandsfontein and Hoedjiespunt) sampled environments that were wetter than LGM–Holocene environments at Nelson Bay Cave (year-round rainfall zone); and (2) considered alongside other evidence from the year-round rainfall zone, wetter conditions across the Pleistocene–Holocene transition at Nelson Bay Cave suggest that climate changes at near-coastal sites may be out of phase with the adjacent interior. Wider application of this approach will provide a semi-quantitative picture of a key paleoenvironmental variable across southwestern South Africa over glacial/interglacial time scales.