Response of coastal California hydroclimate to the Paleocene–Eocene Thermal Maximum

The effects of anthropogenic warming on the hydroclimate of California are becoming more pronounced with the increased frequency of multi-year droughts and flooding. As a past analog for the future, the Paleocene–Eocene Thermal Maximum (PETM) is a unique natural experiment for assessing global and regional hydroclimate sensitivity to greenhouse gas warming. Globally, extensive evidence (i.e., observations and climate models with high pCO₂) demonstrates hydrological intensification with significant variability from region to region (i.e., drier or wetter, greater frequency, and/or intensity of extreme events). Central California (paleolatitude ∼ 42° N), roughly at the boundary between dry subtropical highs and mid-latitude low-pressure systems, would have been particularly susceptible to shifts in atmospheric circulation and precipitation patterns/intensity. Here, we present new observations and climate model output on regional/local hydroclimate responses in central California during the PETM. Our findings, based on multi-proxy evidence within the context of model outputs, suggest a transition to an overall drier climate punctuated by increased precipitation during summer months along central coastal California during the PETM.

1 Introduction

Global warming of a few degrees Celsius over the next century is projected to intensify the hydrological cycle on a range of temporal and spatial scales, manifested primarily by amplified wet–dry cycles (Held and Soden, 2006; Douville et al., 2021). Indeed, just over the last few decades, there has been an increasing frequency in the severity of extremes characterized by compound heat waves and intense drought (Büntgen et al., 2021; Williams et al., 2020; Zscheischler and Lehner, 2022) and/or heavy precipitation and flooding (Liu et al., 2020; Risser and Wehner, 2017). As greenhouse-gas-driven warming continues, such precipitation extremes (wet or dry) are expected to intensify (Stevenson et al., 2022). This is particularly so for California, which receives much of its rainfall from winter systems fueled by atmospheric rivers (ARs), the frequency of which are forecast to decline as the systems shift northward (Simon Wang et al., 2017). The decline in winter precipitation, along with warming, will create more intense droughts even as the potential for extreme precipitation events increases (Vogel et al., 2020; Swain et al., 2018).

Climate model predictions for the intensification of the hydrological cycle are supported by case studies of extreme warming events of the deep past (Carmichael et al., 2017). In particular, the Paleocene–Eocene Thermal Maximum (PETM) has emerged as a unique natural experiment for assessing global and regional hydroclimate sensitivity to greenhouse gas warming (Zachos et al., 2008). Extensive evidence exists for a major mode shift in the local/regional precipitation patterns and intensity (Pagani et al., 2006; Slotnick et al., 2012; Schmitz and Pujalte, 2003; Sluijs and Brinkhuis, 2009; Smith et al., 2007; Handley et al., 2012; Kozdon et al., 2020), including enhanced erosion and extreme flooding in fluvial sections (e.g., Pyrenees and Bighorn Basin) and increased weathering and sediment fluxes to coastal basins (e.g., Bass River, Wilson Lake, mid-Atlantic coast, and Mead Stream, Aotearoa / New Zealand), along with other observations (John et al., 2008; Nicolo et al., 2010; Stassen et al., 2012; Self-Trail et al., 2017; Wing et al., 2005; Kraus and Riggins, 2007; Foreman, 2014).

These observations of regional hydroclimate serve as the basis for climate model experiments forced with proxy-based estimates of ΔpCO₂ for the PETM (i.e., 3×–6× pre-industrial) (Kiehl and Shields, 2013; Carmichael et al., 2016; Zhu et al., 2020). Using such estimates, model simulations show an overall increase in poleward meridional water vapor transport as manifested by a net increase in the evaporation of subtropical regions balanced by higher precipitation of tropical/high latitudes characterizing the “wet-gets-wetter and dry-gets-drier” hydrological response. The latest simulations using high-resolution climate models display several key regional responses, including increased frequency of extreme precipitation events, especially the coastal regions where atmospheric rivers (ARs) are common (Rush et al., 2021). Indeed, observations of high-energy flooding events in SW Europe (i.e., the Pyrenees) during the PETM (Schmitz and Pujalte, 2003) can be explained by increased frequency of North Atlantic ARs contributing landfall in that region. Pacific AR activity as simulated for the PETM also becomes more intense but less frequent along the central California coast by shifting northward with the storm tracks (Shields et al., 2021), not unlike the projections for California in the future (Shields and Kiehl, 2016; Massoud et al., 2019). This pattern is consistent with warming scenarios in general which have weakened zonal wind belts (i.e., the westerlies) that are shifting poleward (Abell et al., 2021; Douville et al., 2021).

Testing the theoretical response of extreme global warming on Pacific ARs and impacts on seasonal precipitation along North America's western coast in general is challenging and still limited by the lack of observations. Here we constrain the regional hydroclimate response along the central California coast during the PETM using several independent proxies (i.e., clay mineralogy, grain size distribution, δ¹³C_org stratigraphy, and leaf wax δ²H_n-alkane isotope records) which are either directly or indirectly sensitive to shifts in precipitation patterns/intensity. These proxies are then compared against Earth system model simulations of the greenhouse-gas-forced changes in regional precipitation (i.e., pattern/intensity). The new records complement data from a previous study (John et al., 2008) and, along with the latest climate modeling experiments, provide a unique case study of the sensitivity of regional hydroclimate to major greenhouse warming.

2 Materials and methods

2.1 Site location

The studied outcrop section, Lodo Gulch, is part of the late Paleocene–early Eocene Lodo Formation located in the Panoche Hills of central California (Fig. 1). During the late Paleocene, the section was situated at a paleolatitude of ∼ 42° N, roughly at the boundary between the dry subtropical highs and mid-latitude low-pressure systems. The Lodo Formation is comprised primarily of siltstone with a relatively low abundance of calcareous microfossils truncated by thin glauconitic sand layers (Brabb, 1983). Depositional facies are consistent with a neritic–bathyal setting along the outer shelf (John et al., 2008).

Figure 1 Paleogeography and location of the Lodo Gulch section (red square) along the central California coast and Bighorn Basin (white square) in the North America continent for reference at 56 Ma. Late Paleocene–early Eocene topography boundary of North America was adapted from Lunt et al. (2017).

2.2 Methods

2.2.1 Bulk organic stable carbon isotopes

Sediment samples used for this study include those originally collected (ca. 28) by John et al. (2008). In addition, new samples (ca. 27) were collected from the upper Paleocene for organic C isotopic analyses (δ¹³C_org) to better establish the pre-PETM baseline. Samples were analyzed in the UCSC (University of California, Santa Cruz) Stable Isotope Laboratory (SiL) using a Carlo Erba Instruments NC 2500 elemental analyzer coupled with a Thermo Scientific Delta Plus XP isotope ratio mass spectrometer (iRMS) via a Thermo Scientific ConFlo III. All samples were calibrated with VPDB (Vienna Pee Dee Belemnite) for δ¹³C and AIR (atmosphere N₂) for δ¹⁵N against an in-house gelatin standard reference material (PUGel). The analytical reproducibility precision is ±0.1 ‰ for δ¹³C and ±0.2 ‰ for δ¹⁵N.

2.3 Grain size analyses

Particle size was measured by laser diffraction using a Beckman Coulter device with polarization intensity differential scatter (PIDS) housed at UCSC (see the Supplement). For each sample, 2 to 5 mg of bulk sediment was powdered and put through a 2 mm sieve following the protocols in Blott et al. (2004). A total of 39 samples was measured, each in duplicate or triplicate to ensure reproducibility.

2.4 Clay assemblages analyses

Sample preparation followed a slightly modified version of Kemp et al. (2016). Roughly 5 to 10 g of sediment was powdered with a mortar and pestle and then placed in a Calgon (sodium hexametaphosphate) solution on a shaker table for 72 h. Samples were sorted through a 63 µm sieve while collecting the fluid with the < 63 µm fraction. The collected fluid and suspended fine fraction (< 63 µm) were allowed to settle for a period determined by Stokes' law to keep < 2 µm size clay particles remaining in suspension. The fluid was then decanted and dried in the oven at 40 °C. Approximately 150 mg clay from each sample were used to prepared oriented mounts for X-ray diffraction (XRD) analysis. A total of 38 clay samples were prepared from the Lodo Formation. The sample residues were measured on a Philips 3040/60 X'Pert Pro X-ray diffraction instrument at UCSC. Clay species (i.e., smectite, illite, kaolinite, and chlorite) were identified based on the peak positions and intensities representing each clay mineral.

2.5 Leaf wax distribution and carbon/hydrogen isotopic composition

Sediment extraction, compound isolation, and compound-specific isotope measurements were conducted following Tipple et al. (2011). Briefly, the sediments were freeze-dried, powdered (∼ 500 g), and extracted with a dichloromethane (DCM) and methanol ratio (2:1; $v/v$) using a Soxhlet extractor. Total lipid extracts were concentrated and then separated by column chromatography using silica gel. N-alkanes were further purified from cyclic and branched alkanes using urea adduction following Wakeham and Pease (2004). N-alkane abundances were determined using a gas chromatograph (GC) with a flame ionization detector (FID). Isotope analyses were then performed using a GC coupled to an iRMS interfaced with a GC-C III combustion system or a high temperature conversion system for δ¹³C and δ²H analyses, respectively. A total of 59 samples was processed with a fused silica, DB-5 phase column (30 m × 0.25 mm internal diameter, 0.25 µm film thickness) with helium as the carrier at a flow of 1.5 mL min⁻¹. The GC oven temperature program was 60–320 °C at 5 °C min⁻¹ and isothermal for 30 min. A Thermo Trace mass spectrometer (MS) was used for detection with the mass spectrometry scanning from 50–800 $m/z$ or exclusively $m/z$ of 191, 217, 218, 370, 372, 386, and 400 for single-ion monitoring. Biomarkers were identified by the elution time and mass spectra of in-house petroleum standards with published biomarker distributions (Peters et al., 2005).

δ¹³C and δ²H values are expressed relative to Vienna Pee Dee Belemnite (VPDB) and Vienna Standard Mean Ocean Water (VSMOW). Individual n-alkane isotope ratios were corrected to n-alkane reference materials (for δ¹³C, C₂₀, C₂₅, C₂₇, C₃₀, and C₃₈ of known isotopic ratio and for δ²H; “Mix A” from Arndt Schimmelmann, Indiana University) analyzed daily at several concentrations. In addition, H₂ reference gas of known isotopic composition was pulsed between sample n-alkane peaks to confirm if the normalizations were appropriate. Standard deviations (SDs) of n-alkane reference materials were ±0.6 ‰ for δ¹³C and ±6 ‰ for δ²H.

2.6 Earth system models

Climate simulations from two models were used in this paper for (1) comparison with leaf wax proxy data and (2) extreme events analyses. (1) The water-isotope-enabled Community Earth System Model version 1.2 (iCESM1.2) simulates changes in climate and water isotopic composition during the PETM (Zhu et al., 2020) with a horizontal resolution of 1.9 × 2.5° in atmosphere and land and a nominal 1° in the ocean and sea ice components. Water isotope capabilities have been incorporated into all the components of the Community Earth System Model (CESM) 1.2 (Brady et al., 2019) which include the Community Atmosphere Model, version 5 (CAM5), for the atmosphere; the Parallel Ocean Program, version 2 (POP2), for the ocean; the Community Land Model, version 4 (CLM4), for the land; the River Transport Model (RTM) for river flow; and Community Ice Code, version 4, for sea ice. All simulations were run with the identical boundary conditions (including early Eocene paleogeography, land–sea mask, vegetation distribution, and pre-industrial (PI) non-CO₂ greenhouse gas concentrations, soil properties, natural aerosol emissions, solar constant, and orbital parameters), following the DeepMIP protocol (Lunt et al., 2017), and differ only in atmospheric CO₂ concentration. Crucially, the models with reduced latitudinal temperature gradients (e.g., Geophysical Fluid Dynamics Laboratory (GFDL) and CESM) more closely reproduce proxy-derived precipitation estimates and other key climate metrics (Cramwinckel et al., 2023). Increased climate sensitivity with warming and cloud feedback in CESM1.2 over earlier models improved water vapor sensitivity. (2) Using the same CESM1.2 framework, high-resolution (0.25°) simulations were conducted with forced sea surface temperatures (SSTs) and active atmosphere and land components (CAM5 and CLM4). RTM was run at 1° resolution, and forced SSTs were calculated from consistent 2° fully coupled PETM simulations (see details in Rush et al., 2021, and references therein). The much higher horizontal resolution in the atmosphere enables improved the simulation of the extreme events. Hourly, daily (CAM5), and monthly (iCESM1.2) temporal resolution precipitation outputs from both sets of climate simulations were utilized in this paper, with 100 years taken from the equilibrated iCESM1.2 simulations and 15 years from the forced SST high-resolution CAM5 simulations.

3 Results

3.1 Bulk organic and n-alkane stable carbon isotopes

A carbon isotope excursion is present in both bulk organic (Fig. 2a) and carbonate-based δ¹³C records (John et al., 2008) across the P−E boundary, marking the PETM onset of the Lodo section. The terrestrial leaf wax n-alkane records all capture the carbon isotope excursion (CIE) with a pattern that roughly parallels the other published records (i.e., planktonic foraminifera) (John et al., 2008), though is much less noisy than the bulk δ¹³C_org record, which is not unexpected given the potentially variable composition of the bulk organic matter. The magnitude of the Δδ¹³C_n-alkane is roughly −4 ‰ (average of n-C₂₇, n-C₂₉, and n-C₃₁) at the onset of the CIE, followed by a gradual recovery that is truncated at the disconformity between 20.3 and 23.5 m (coincides with nannofossil biozone boundary NP10 and NP11), thus marking the top of the PETM body (Fig. 2a). The disconformity coincides with a global sea level regression (John et al., 2008). Following the recovery, above the disconformity, the mean δ¹³C_n-alkane is depleted relative to the pre-PETM baseline, as observed in some other PETM sections (Cui et al., 2021; Garel et al., 2013; Hasegawa et al., 2006).

Figure 2 Terrestrial higher plant leaf wax n-alkane δ¹³C and δ²H records for the Lodo Gulch section, central California. The shaded area represents the bounds of the CIE/PETM (a) bulk organic carbon isotope record. (b, c) Leaf-wax-compound-specific carbon/hydrogen isotope records in n-C₂₇ (yellow square), n-C₂₉ (red closed circle), n-C₃₁ (blue triangle), and (d) n-alkane carbon preference indices (CPIs).

3.2 Hydrogen isotopes

The leaf wax δ²H_n-alkane values range from −150 to −213 ‰ over the entire sampled section that has an initial decrease of 25 ‰ (from −150 to −175 ‰ in C₂₉) just prior to the CIE onset and then followed by a slight rise (∼ 6 ‰) right after the onset. The relatively invariable δ²H_n-alkane through the PETM is punctuated with two brief intervals of more negative values (−202 ‰ at 6.26 m and −213 ‰ at 22 m). The second larger anomaly coincides with the disconformity (related to local sea level regression). The post-PETM δ²H_n-alkane values are on average lower than for the upper Paleocene/PETM, although analytical errors may bias the values. Given the limited number of samples to establish a baseline for the upper Paleocene, the significance of the pre-CIE (and post-CIE) shifts/anomalies in δ²H_n-alkane should be considered with some caution. Several other sections do show pre-CIE shifts, both positive and negative, and typically an enrichment with the CIE (Handley et al., 2008, 2011; Jaramillo et al., 2010; Tipple et al., 2011). Such minor changes likely reflect unconstrained orbital influences on regional precipitation (Rush et al., 2021; Campbell et al., 2024), especially considering the variable direction of change from location to location.

3.3 Clay assemblage and grain size

Clay assemblages and particle grain size should to some extent be influenced by regional hydroclimate. At Lodo, the clay assemblages (Fig. 3) are dominated by smectite throughout. The minor clay components of illite and chlorite show several spikes relative to smectite within the lower (8 to 10 m) and upper CIE (∼ 19 m), whereas the ratio of kaolinite gradually increases (0.5 to 1.5) only over the upper portion of the CIE (10 to 20 m). A delayed rise in kaolinite has also been observed in a few other PETM sections, whereas some show an immediate rise (Tateo, 2020; Gibson et al., 2000). The smectite concentration and $\mathrm{kaolinite}/\mathrm{smectite}$ ratio remain high in the post-PETM interval. The late Paleocene, with a relatively coarse sand size, shows slight spikes of kaolinite associated with other minerals. Grain size, largely for silt and clay, shows a distinct shift toward finer fractions (i.e., clay) with the onset of the CIE (Fig. 3e).

Figure 3 Integrated C isotope and clay assemblage records for the Lodo Gulch section, central California. (a) Bulk organic carbon isotope. (b, c, d, e) Clay assemblage ratios. (e) The decreasing mean particle size (D50: 50 % of the total particle size in sediments) corresponds to CIE onset.

3.4 Earth system model simulations

We obtained and processed temperature and precipitation outputs from two community Earth system/climate models: the isotope-enabled iCESM1.2 and high-resolution CAM5 models (daily precipitation over 15 years) forced by a range of greenhouse conditions (1×, 3×, 6×, and 9× pCO₂ pre-industrial), both with Eocene paleogeography. For comparisons with observations, we used the output from the 3× to 6× pCO₂ simulations which best replicated the observed SSTs (ΔSST) for the pre-PETM and PETM (Zhu et al., 2020). Overall, monthly winter precipitation for the study region decreases (∼ 30 %) during the PETM in both simulations, with a slight increase in the summer (Figs. 4, 5). CAM5 output shows a modest decrease in the mean annual precipitation, with significant seasonal shifts during the PETM (Fig. 5a). Seasonal changes in the monthly averaged δ¹⁸O and δ²H from the mean monthly precipitation in iCESM1.2 of central California are consistent with CAM5. On average, the δ²H_precip increases by ca. 5 ‰–10 ‰ from pre-PETM to PETM, especially in the winter/spring, with a smaller shift in summer/fall (1 ‰–2 ‰) (Fig. 4a, b, c). The extreme value index (ξ), a representation of the distribution of the exceedance right tail (see the Supplement), shows a small but statistically robust increase in the wet extremes of winter (December–February, DJF), with a significant increase in summer (June–August, JJA) during the PETM (Fig. 5b).

Figure 4 Seasonal and monthly meteoric precipitation amounts (mm d⁻¹) and H-isotopic composition for the North Pacific/western North America as simulated with the iCESM1.2 (Zhu et al., 2020) under pre-PETM (3× in black) and PETM (6× in red) pCO₂ forcing. Panels (a) and (b) show Δδ²H_precip between pre-PETM (3×) and PETM (6×) in winter (DJF) and summer (JJA). Panel (c) shows the annual seasonal cycle of δ²H_precip for central California (pre-PETM in black; PETM in red). Mean daily precipitation rate difference for (d) winter and (e) summer between pre-PETM (3×) and PETM (6×). Panel (f) shows the annual seasonal cycle of daily precipitation rate for central California. Values represent the area-weighted average over 4° × 4° box bounding the study site.

Figure 5 High-resolution CAM5 model output (Shields et al., 2021) of (a) mean monthly precipitation for central California over 15 model years for the late Paleocene (LP) under low pCO₂ (680 ppmv) and for the PETM under high pCO₂ (1590 ppmv). (b) Extreme value index (ξ) comparison of mean monthly precipitation in winter and summer of the central California region based on the CAM5 precipitation outputs.

4 Discussion

4.1 Hydroclimate response from model simulations

In all model simulations of the PETM forced with higher pCO₂ (e.g., 3× to 6× pre-industrial), the hydrological cycle intensifies, as manifested by increases in global mean precipitation and meridional vapor transport (Kiehl and Shields, 2013; Kiehl et al., 2018; Carmichael et al., 2016; Zhu et al., 2020). Regionally however, the magnitude and even the sign of precipitation change can differ considerably from global means (Carmichael et al., 2016, 2017). This is most evident in the latest low- and high-resolution model simulations of the PETM (Zhu et al., 2020; Shields et al., 2021). For central California, the simulations yield an overall decline in the mean annual precipitation that is mainly due to a notable decline in winter precipitation with only a slight increase in summer (Fig. 4). This pattern is produced by both the water-isotope-enabled iCESM1.2 and the higher-resolution CAM5 with an overall shift to lower-amplitude seasonal cycles (i.e., drier winter/spring and a slightly wetter summer) (Figs. 4, 5). This seasonal wet–dry shift appears to be driven in part by a pronounced northward shift in atmospheric rivers (ARs) in winter along the North American Pacific coastline (Shields et al., 2021). As ARs deliver most of the winter precipitation to the mid-latitude Pacific coast, less frequent occurrences result in drier winters during the PETM. Moreover, the extreme value index (ξ) shows a small but statistically robust increase in winter (DJF) wet extremes, with a significant increase in the probability of summer (JJA) wet exceedance during PETM (Fig. 5b). The latter might be due to elevated tropical storm activity along the Pacific coast during PETM (Fig. S7; Kiehl et al., 2021).

4.2 Hydroclimate response from observations

4.2.1 Sedimentation rate, clay assemblages, and grain size distribution

Arguably, the collection of observations from Lodo (i.e., sediment flux, clay assemblages, and leaf wax δ²H isotopes), within limitations, appears to be mostly consistent with the model output. Starting with siliciclastic sedimentation, rates should be highly susceptible to a major shift in hydrologic conditions as changes in the seasonality of precipitation (along with vegetative cover) would impact rates of erosion and sediment transport. The coarse resolution of stratigraphic control at Lodo does limit the ability to constrain changes in sedimentation rates in detail. However, just considering the thickness of the CIE (∼ 10 m), one could argue for a shift toward higher seasonality of precipitation with overall drier conditions, as suggested by John et al. (2008).

The other constraints on regional precipitation would also support a shift toward drier conditions. For example, an increase in the relative abundance of kaolinite fluxes has been widely observed across the CIE onset in many PETM sections from mid- to high latitudes (Tateo, 2020; Gibson et al., 2000) and interpreted as evidence of a major mode shift in local hydroclimates. In contrast, the clay mineralogy (Fig. 3) for the Lodo Formation is dominated mainly by smectite at the onset of the PETM, consistent with seasonal wet/dry cycles under warm conditions (Gibson et al., 2000). A subtle increase in the $\mathrm{kaolinite}/\mathrm{smectite}$ could be interpreted as evidence of higher humidity (Gibson et al., 2000). However, the skewed grain size distribution of clay sediments at around 8 m, coinciding with illite/smectite peaks (Figs. 3e and S1), indicate higher fluvial velocity and increased erosion, as observed elsewhere (Chen et al., 2018; Foreman et al., 2012; Foreman, 2014). For example, along the mid-Atlantic margin it appears that the kaolinite might have been exhumed from local Cretaceous laterites (Lyons et al., 2019). The regionally enhanced physical weathering and erosion could be related to an increase in the frequency of episodic wet/dry extremes during the PETM, as simulated by models.

4.2.2 Precipitation and leaf wax δ²H_n-alkane

The Lodo leaf wax δ²H_n-alkane record is, at first glance, somewhat equivocal in terms of the response of local hydroclimate. In theory, terrestrial higher plant δ²H_n-alkane should provide insight into the changes in regional precipitation amounts/source, particularly for major mode shifts (Jaramillo et al., 2010; Pagani et al., 2006; Tipple et al., 2011). In some PETM records, δ²H_n-alkane significantly increases, consistent with the effects of higher T on water isotope fractionation. For example, in the Arctic, δ²H_n-alkane records show a positive excursion of 55 ‰ at CIE onset, consistent with higher T, a reduced meridional temperature gradient, and a decrease in isotope distillation during vapor transport (Pagani et al., 2006). However, there are notable exceptions. In some subtropical/mid-latitude sites, δ²H_n-alkane decreases (ca. ∼ 20 ‰) across the onset of PETM (Handley et al., 2008, 2011; Jaramillo et al., 2010; Tipple et al., 2011). In comparison, Lodo δ²H_n-alkane displays a comparatively muted response, showing a slight ²H enrichment in the main body PETM, followed by several anomalous shifts toward more negative values (Fig. 2c).

Given the robust evidence for mode shifts in hydroclimate elsewhere during the PETM, does the relatively stable Lodo δ²H_n-alkane record necessarily support a local/regional stable hydroclimate (i.e., in conflict with the modeling and other observations)? As H-isotope fractionation in plants is related to photosynthetic pathways, source water availability, and atmospheric humidity (Sachse et al., 2012; Tipple et al., 2015), it is possible that local shifts in meteoric water isotope composition were offset by another influencing factor(s). Regarding photosynthetic pathways, along the west coast of North America, no detailed records of vegetation response have been generated for the PETM. Still, for the late Paleocene and early Eocene intervals, Korasidis et al. (2022) found little deviation in the Köppen–Geiger climate type (i.e., Mediterranean) within the central California region. This evidence, along with the lack of change in the n-alkane average chain length (ACL) in the Lodo section (Fig. S5), would suggest no major changes in vegetation assemblages during the PETM. Another factor, a reworking of the Paleocene terrestrial organic matter (e.g., Tipple et al., 2011), could possibly dampen of isotopic n-alkane signals at Lodo, although the CPIs and the leaf wax carbon isotopes would suggest minimal reworking of the n-alkanes, as opposed to other coastal PETM sites where the evidence for reworking is robust (e.g., Lyons et al., 2019; Hollingsworth et al., 2024). As such, if we assume the δ²H_n-alkane record reflects only on changes in local meteoric waters, then the observed modest change in δ²H_n-alkane values at Lodo could be interpreted in several ways in terms of T-related changes on isotope fractionation that were offset by changes in the dominant season of precipitation and/or vapor sources and distance of transport. For example, a shift in precipitation between winter and late summer/fall could offset the effects of warming, assuming a shift from a proximal (north or central Pacific) to a more distal (Gulf of Mexico) source of vapor (Hu and Dominguez, 2015). At ground level, stronger evapotranspiration during biosynthesis can isotopically be offset by external water source availability (i.e., seasonal precipitation). Local/regional groundwater table variations caused by a hydrological change would also affect the source water use efficiency of plants since surface water tends to be more depleted in some perennial species after intense storms in the groundwater (Hou et al., 2008; Krishnan et al., 2014). Hydrogen isotope fractionation in plants can also be biased by seasonal shift in regional vegetation growth regime. For example, leaf wax lipids from terrestrial plants usually record hydrological conditions earlier in the season rather than fully integrating the entire growing season (Hou et al., 2008; Tipple et al., 2013). Finally, episodic extremes in precipitation may dominate the hydrogen isotopic composition of the leaf wax (Krishnan et al., 2014). If soil water is derived mainly from extreme events during the growth season, then the lack of a major shift in the Lodo δ²H_n-alkane record with the onset of the PETM could reflect a combination of more ²H-depleted precipitation delivered by seasonal storms offset by warming-induced ²H-enrichment in leaf water.

4.3 Comparison of leaf water and modeled δ²H

Assuming leaf wax δ²H is primarily influenced by local meteoric water, how does the seasonal distribution of precipitation influence the bulk signal? In iCESM1.2 simulations with increasing pCO₂ (i.e., 3× to 6× pre-industrial) and SSTs, the seasonal shifts in δ²H of mean monthly precipitation from pre-PETM to PETM are significant. During the winter, as precipitation amounts decline, δ²H_precip increases by 10 ‰ while decreasing by ∼ 1 ‰ to 5 ‰ during late summer/fall in central California (Fig. 4). To estimate how this seasonal change in δ²H_precip and precipitation amount influences leaf water δ²H, we applied a leaf wax proxy model (see the Supplement) which computes the combined effects of changes in seasonal precipitation and growing season length. The model shows leaf water δ²H is enriched ca. 4 ‰ to 7 ‰ from pre-PETM to PETM. Arguably, this would be consistent with the minor ∼ 6 ‰ enrichment observed in the Lodo record at the onset of the PETM. We also examined other sites for comparison of the predicted leaf water differences from pre-PETM to PETM with the fossil leaf wax Δδ²H (see the Supplement) and find a similar pattern in other mid-latitude sites. Other factors to consider include precipitation source waters and a shift in a mixing endmember between proximal and distal sources of water in the coast. For example, with a summer shift in the source water from the Pacific to subtropics (i.e., summer monsoons), the effect of increasing distance and distillation would isotopically deplete vapor (Hu and Dominguez, 2015), thus offsetting the temperature-related enrichment of local δ²H_precip. In addition, infrequent but high-intensity tropical cyclones during the PETM (Kiehl et al., 2021) would tend to deliver relatively depleted precipitation (i.e., a more negative δ²H) during summer months.

Finally, a related record that might indirectly reflect on precipitation amount (i.e., atmospheric humidity) is the magnitude of the CIE recorded by leaf wax δ¹³C_n-alkane. Recalcitrant higher plant leaf wax n-alkane carbon isotope ratios (n>25 with an odd-over-even preference) mainly reflect the carbon source (Diefendorf et al., 2010). However, photosynthetic carbon isotope fractionation (Δ_p) is sensitive to atmospheric pCO₂ variations, generally increasing with rising concentrations assuming a constant photosynthetic fractionation factor and humidity (Diefendorf et al., 2010). The δ¹³C_n-alkane of Lodo section displays a sharp negative shift of ca. 4 ‰ (average of n-C₂₇, n-C₂₉, and n-C₃₁) across the onset of CIE (Fig. 2b), which is consistent with global mean atmospheric CIE (Sluijs and Dickens, 2012) but generally smaller than observed in other leaf wax records (Handley et al., 2008, 2012; Jaramillo et al., 2010; Pagani et al., 2006; Tipple et al., 2011). The smaller δ¹³C_n-alkane CIE recorded in Lodo would be consistent with a reduction in local humidity which should reduce the magnitude of Δ_p during photosynthetic carbon fixation.

5 Summary and conclusions

Many sections globally exhibit evidence, often striking, of significant shifts in local hydroclimate at the onset of the PETM consistent with model simulations (e.g., 2× CO₂). These same models also simulate an overall decrease in winter precipitation for the central California coast due in large part to a reduction in AR frequency (Shields et al., 2021). While not as striking, the collection of observations from the central California Lodo Gulch section would support a modest reduction in precipitation (i.e., MAP, mean annual precipitation) during the PETM, along with the possibility of an increase in the frequency of extreme precipitation events. This transition toward greater aridity and precipitation extremes is not unlike the forecasts for much of California over the coming centuries due to anthropogenic warming.