Chert oxygen isotope ratios are driven by Earth's thermal evolution

To explore the dependence of δ¹⁸O_chert on silica transformation kinetics, we investigated the covariation between authigenic quartz δ¹⁸O (δ¹⁸O_a-qtz) and f_sil in cherts and silicified shales from the 550- to 525-Ma-old (39) Lijiatuo section in South China. The bulk δ¹⁸O values vary between 12.6 and 24.0‰, and f_sil ranges from 1 to 70 wt% across the ca. 140-m stratigraphic section. Calculated median δ¹⁸O_a-qtz values (derived via mass balance; SI Appendix, Supporting Information 1) range from 12.6 to 28.7‰ and covary with f_sil (trend in Fig. 2A; r² = 0.56 and P < 0.01). Secondary Ion Mass Spectrometry (SIMS) analyses of microcrystalline quartz confirm their high δ¹⁸O in detrital-rich samples; low and intermediate δ¹⁸O values represent pure or partially overgrown detrital quartz grains (SI Appendix, Supporting Information 1 and Fig. 2A). Using the temperature–fractionation calibration of Sharp et al. (1) and assuming porewater δ¹⁸O = −1.1‰ yields temperatures of 58 °C to 189 °C. Both the absolute values and the range of temperatures are implausible for seawater. Alteration with meteoric water or reequilibration with fluids at elevated temperatures could shift δ¹⁸O to lower values but cannot explain the dependence of δ¹⁸O_a-qtz on f_sil . Instead, the data are best explained by a change in the rates of silica polymorph transformation during diagenesis, such that quartz precipitates earlier at higher f_sil . This follows the long-standing interpretation of silica polymorph distributions in the Californian Monterey Formation, where the opal-CT/quartz transformation occurs at shallower burial depth and lower temperatures as f_sil increases (34). The underlying mechanism is thought to be reduced silica solubility at high f_sil that causes higher initial crystallinity in the intermediate opal-CT, such that the time required for crystallographic ordering before transformation to quartz is reduced (34). In summary, the covariation of δ¹⁸O_a-qtz with f_sil confirms that δ¹⁸O_chert is predominantly controlled by the suite of parameters that determine the kinetics of silica diagenesis and thereby the temperature and timing of the final silica phase transformation.

Lijiatuo section δ¹⁸O_a-qtz-inferred temperatures do not correlate with age (as proxied by height in the profile) and are much lower than peak diagenetic temperatures, estimated by Raman spectra of carbonaceous material to ~280 °C (Fig. 2B; Materials and Methods). Our data thus extend, by >500 Ma, the observation from the Miocene Monterey Formation sediments that δ¹⁸O_a-qtz remains unchanged following quartz formation (30). The preservation of the δ¹⁸O_a-qtz–f_sil relationship testifies to a pristine δ¹⁸O_chert signature set during the final opal-CT/quartz transformation—despite later higher-grade diagenesis and the potential for modern alteration. This observation emphasizes the robustness of δ¹⁸O_chert as a recorder of diagenetic conditions in sedimentary cherts and further stresses the importance of transformation kinetics and the prograde T–t pathway in setting δ¹⁸O_chert.

We have developed a 1D reaction–advection–diffusion model to describe the diagenesis of silica polymorphs (full details in Materials and Methods and SI Appendix, Supporting Information 2). At the heart of the model is a Si mass balance, which couples the dissolution and precipitation of three silica polymorphs (opal-A, opal-CT, and quartz) through a porewater nexus. In the model, opal-A is deposited onto an accumulating sediment column and progressively transforms to opal-CT and quartz as depth and temperature increase. The rates of SiO₂ polymorph dissolution and reprecipitation are governed by kinetic rate constants and the degree of local porewater under- or oversaturation with respect to each phase. Temperature sets the magnitude of oxygen isotope fractionation and influences diagenetic transformation rates via Arrhenius formulations of the kinetic constants of dissolution and reprecipitation. Temperature also influences the solubilities of the three SiO₂ phases, the diffusion coefficients for porewater dissolved Si, and thus the thermodynamic driving force for any given reaction. The model includes the acceleration of the opal-CT to quartz transition with increasing detrital mineral fraction, which is empirically calibrated against the Lijiatuo section data (SI Appendix, Supporting Information 2), and successfully reproduces the quartz precipitation depth and δ¹⁸O_chert values in Sea of Japan sediment cores (27).

Our model demonstrates that the fundamental control on the offset between δ¹⁸O_chert and δ¹⁸O_SW is the temperature window over which quartz precipitates and that this temperature is partially decoupled from seawater temperature. There is not one single temperature or depth at which opal-CT transforms to quartz (32). Rather, it depends on the suite of parameters that control the rate of progression along the entire SiO₂ diagenetic pathway that variably influence δ¹⁸O_chert (Fig. 3). The model reveals that ambient pore fluid δ¹⁸O is largely fluid-buffered, and hence, variation in water–rock ratios is not an important driver of δ¹⁸O_chert. Remarkably, heat flow (Q, W m⁻²) has substantial leverage on δ¹⁸O_chert (Fig. 3B), which means that δ¹⁸O_chert paleothermometry must account for paleo heat flow. Model sensitivity analyses indicate that besides heat flow, δ¹⁸O_SW, and the decoupled T_SW effect, other parameters have little leverage over the secular trend in δ¹⁸O_chert (SI Appendix, Supporting Information 3). Below, we show that the heat flow effect reconciles the early Paleozoic δ¹⁸O_chert with only moderately depleted δ¹⁸O_SW and moderate seawater temperatures and that low Archean δ¹⁸O_chert requires neither high ocean temperatures nor very low δ¹⁸O_SW, given reasonable estimates of paleo heat flow.

Heat flow through ocean crust in the Phanerozoic can be quantified due to the well-established decrease in heat flow with crustal age (42) and the availability of reconstructions of crustal age distributions through time (43). In contrast with Precambrian time, most Phanerozoic chert forms in deep ocean settings. The average ocean crustal age in the early Paleozoic was ca. 20 Myr, increasing to ca. 50 Myr for the late Cenozoic (43). Applying the age–heat flow relationship (42) yields Paleozoic average heat flow values of ca. 105 mW m⁻² and Cenozoic values of 80 mW m⁻² (42). All else being equal, the heat flow effect from different ages of ocean crust alone accounts for ca. 2‰ of the 8‰ lower δ¹⁸O_chert of Cambrian age relative to chert from the Eocene (Fig. 1). Considering heat flow thus solves the dilemma of δ¹⁸O_chert-inferred seawater temperatures above the metazoan lethal threshold without necessitating strongly ¹⁸O-depleted oceans (44). Nonetheless, the lack of high δ¹⁸O_chert values during the well-documented Paleozoic icehouses, combined with our demonstration of the insensitivity of δ¹⁸O_chert to postprecipitation alteration, suggests that early Phanerozoic oceans were indeed moderately ¹⁸O depleted.

Over longer timescales, a reduced sensitivity of δ¹⁸O_chert to changes in T_SW precludes the hot Archean ocean hypothesis as the sole explanation for low δ¹⁸O_chert. A temperature reduction (ΔT_SW) from 70 to 5 °C at Archean heat flow yields an increase in δ¹⁸O_chert of 9‰ (Fig. 3C) compared with 16‰ from direct application of the temperature–fractionation (ε^18/16O) relationship for the silica–water system. This behavior results from two mechanisms. First, the nonlinear relationship between ε^18/16O and temperature (1) (SI Appendix, Fig. S7) means the same ΔT_SW translates to a smaller isotope difference at elevated diagenetic temperatures. Second, the offset between diagenetic temperatures and T_SW is reduced at higher T_SW (Fig. 3A), narrowing any given ΔT_SW.

Plausible decreases in heat flow explain a large fraction of the 15‰ δ¹⁸O_chert increase across Earth history. Although detailed tectonic reconstructions are lacking for most of the Precambrian, diverse lines of evidence indicate that heat flow through Archean crust was around 2× higher (45–47) than the modern value of ca. 80 mW m⁻² (42, 45, 48, 49). A decrease of heat flow from 160 to 80 mW m⁻² corresponds to an increase in δ¹⁸O_chert of ca. 5.4‰ at T_SW = 5 °C, of 3.3‰ at T_SW = 30 °C, and of 1.3‰ at T_SW = 70 °C (Fig. 3 A and C). This control emerges from the kinetics of silica diagenesis: at lower heat flow, the opal-CT/quartz transformation happens later in time and at greater depth but at a lower temperature and thus with a larger silica–water fractionation (Fig. 3A). The heat flow effect on δ¹⁸O_chert was previously overlooked and removes the requirement for extreme changes in δ¹⁸O_sw to explain ¹⁸O-depleted cherts, pointing the way toward a middle-ground scenario that satisfies all biological, geochemical, and paleoclimatological constraints without invoking pervasive alteration.

A 2× decrease in heat flow, in combination with either no temperature change or temperature reductions of 25 °C or 65 °C, produces shifts in δ¹⁸O_chert of ca. 5, 10, or 14‰, respectively (Fig. 3D). Since we discount metamorphic or meteoric alteration, δ¹⁸O_SW is the sole remaining explanation for the residual shift in the δ¹⁸O_chert record (−10, −5, and −1‰, respectively). Thus, the scenario of hot Archean oceans with near-modern δ¹⁸O_sw (50) can be mathematically salvaged if considered in conjunction with a ≥2× decrease in heat flow, but it remains biologically implausible and inconsistent with Precambrian glaciations. Beyond this, it is also incompatible with relatively constant δ¹⁸O_chert ranges (cf. box-and-whisker plots in Fig. 1). Our model predicts that for a ΔT_SW of 65 °C (from 70 °C to 5 °C), the δ¹⁸O_chert range should substantially increase (Fig. 3C). This behavior results from silica diagenesis occurring rapidly, at shallow depths, in all high-T_SW scenarios regardless of local heat flow, depositional setting, or sediment properties. Conversely, at lower T_SW, the leverage of local variations in these properties on the overall silica transformation rate is increased, producing a broader distribution of δ¹⁸O_chert. We thus infer that the consistency in the range of δ¹⁸O_chert through time is evidence for consistency in T_SW.

The seawater ¹⁸O budget is itself a function of solid Earth dynamics (51, 52) and reflects the global ratio of high-T to low-T alteration processes. High-temperature hydrothermal fluids are the dominant source for ¹⁸O because equilibrated fluids are close to a basalt δ¹⁸O value of ca. 5.8‰. In contrast, low-T alteration products (e.g., clays) preferentially remove ¹⁸O from the hydrosphere. The connection between the seawater ¹⁸O budget and solid Earth dynamics thus arises by changing the high-T to low-T alteration ratio (14, 51). This can be achieved by altering parameters including crustal spreading rates or CO₂ degassing rates, the rate of continental emergence, plate tectonic regime, or global sea level. All these factors are known to influence the rate and temperature of oxygen-exchanging water–rock interactions (12, 14, 51). A recent numerical model that explicitly links the fluxes of ¹⁸O to and from the hydrosphere with the thermal evolution of the solid Earth (52) suggests seawater δ¹⁸O at 3.0 Ga of ca. 5 to 10‰ lower than today. These model predictions combine with the 1.3 to 5.4‰ decrease in δ¹⁸O_chert resulting from a 2× higher Archean heat flow (see above) to potentially explain essentially the entirety of the chert oxygen isotope record (Fig. 1). At a minimum, the thermal effects on δ¹⁸O_chert sum to at least ca. 6.3‰, such that Archean seawater temperatures >30 °C are precluded (Fig. 3D). Overall, we conceptualize that δ¹⁸O_chert encodes a dampened T_SW signal and both direct (heat flow) and indirect (δ¹⁸O_sw) effects of Earth’s thermal evolution.

There are also controls over the ¹⁸O fluxes that do not directly reflect solid Earth dynamics. More efficient reverse weathering reactions—the result of a high Precambrian ocean Si inventory (53)—would act to increase ¹⁸O consumption directly by low-temperature mineral precipitation but also indirectly by increasing the net rate of continental silicate weathering—itself an ¹⁸O consuming reaction—that has declined through time according to geochemical models (53, 54). The continental surface has become gradually more covered by higher δ¹⁸O sedimentary lithologies (55) with less leverage over the ocean ¹⁸O seawater budget. Declining rates of reverse weathering and increasing continental sedimentary coverage thus both act to raise δ¹⁸O_SW through time. Overall, we infer that a changing balance between high-T and low-T water–rock interactions through time has increased δ¹⁸O_SW by up to 9‰ with the precise number depending on the Archean climate scenario. By accounting for heat flow, we are able to identify a compromise solution to the long-standing T_SW/δ¹⁸O_SW debate, without having to invoke ubiquitous alteration or fortuitously monotonic climate evolution.

We show that the δ¹⁸O_chert record is a largely pristine record of diagenetic conditions and is sensitive to heat flow. This finding implies that the long-term, near-monotonic increase in δ¹⁸O_chert is an expression of Earth’s thermal evolution. Accounting for heat flow helps explain the ¹⁸O-depleted early Phanerozoic and Archean record without needing to invoke implausibly hot seawater or radically different δ¹⁸O_sw. We show that ca. 5‰ of the 15‰ Archean-Cenozoic difference in δ¹⁸O_chert is directly attributable to a change in diagenetic temperatures caused by decreasing basal heat flow and in situ crustal heat production (Fig. 3). This effect combined with moderately depleted δ¹⁸O_sw of ca. −5 to −9‰—in line with a recent δ¹⁸O_sw evolution model explicitly forced by reconstructions of Earth’s crustal and thermal evolution (52)—reconciles the chert δ¹⁸O record with geological evidence for “icehouse” periods and high CO₂ sequestration fluxes in the late Archean (22). Geographic variation in seawater temperatures, chert depositional environments, sedimentation rates, and variations in detrital mineral contents contribute to the δ¹⁸O_chert scatter at any given time. Seawater temperatures leave only a dampened signal in the record of δ¹⁸O_chert.

The framework we advance here is applicable to the oxygen isotope records of carbonates and other chemical sediments that recrystallize (i.e., undergo dissolution–reprecipitation reactions) or otherwise exchange with pore fluids during diagenesis. For example, foraminiferal tests have been shown to equilibrate their δ¹⁸O with pore fluids at elevated temperatures during burial (56, 57). Since the extent of recrystallization and exchange are dependent on temperature and time, the δ¹⁸O of other chemical sediments will also respond to changes in heat flow, although they will need specifically tailored assessments of their diagenetic pathways to fully exploit. Broadly, the similarity of δ¹⁸O among carbonates, cherts, and other archives through geological time (13, 14) implies similar sensitivities to declining heat flow. Our framework also hints at future directions for empirical paleoclimate reconstructions when constraints on paleo heat flow and δ¹⁸O_sw are available. Overall, the records of δ¹⁸O from chemical sediments are strongly suggestive of a system driven by forcings that vary on billion-year timescales. Declining heat flow and increasing δ¹⁸O_sw as the Earth cools provide such a driver in a way that the more capricious climate system cannot.

We analyzed cherts and silicified shales from the “Lijiatuo” section (28°24.519 N and 110°27.926 E) on the Yangtze Platform, China. We determined major element concentrations on fused bulk rock samples by X-ray fluorescence (XRF) with a Philips Panalytical PW2400. Sample X-ray diffraction (XRD) patterns were generated with a Panalytical Empyrean, and mineral phase composition and quantification were evaluated using Autoquan. The total organic carbon content (TOC) was determined on decarbonated samples (20% HCl) using an NA1500 elemental analyzer. Oxygen isotope ratios of bulk rocks were determined by laser fluorination gas source mass spectrometry (58), in which oxygen is liberated from samples in the presence of 10 mbar F₂ with a SYNRAD 50 W CO₂ laser. O₂ is purified from other reaction products in a series of liquid nitrogen cold traps and molecular sieves and ultimately introduced in continuous flow mode into a Thermo MAT 253 mass spectrometer. Oxygen isotopes were simultaneously measured on m/z = 32 (¹⁶O¹⁶O) and m/z = 34 (¹⁶O¹⁸O) for samples and the NIST NBS 28 standard and are reported as δ values on the VSMOW scale. In situ analyses of quartz grains was performed in epoxy-embedded samples using a Cameca 1280-HR SIMS with a primary ¹³³Cs⁺ beam to sputter O-ions. A low-energy electron gun was used to minimize surface charge accumulation. NIST NBS 28 (δ¹⁸O_VSMOW = +9.57±0.10] was used to correct for instrumental mass bias and anchor the oxygen isotope ratios to the VSMOW scale. Results have an estimated reproducibility of ±0.35‰. Raman spectra of carbonaceous material were acquired using a Horiba Jobin Yvon HR800-UV spectrometer connected to an Olympus BX41 microscope operating in 180° backscattering geometry using a 488-nm Ar⁺ laser for excitation. On-sample laser power was attenuated by density filters to <0.5 mW. Calibration of Raman bands was achieved using the 520.4 cm⁻¹ silicon band, and spectra converted to an estimate of peak diagenetic temperature based on a temperature-dependent change in position and width of Raman bands in the first-order region of carbonaceous matter Raman spectra.

δ¹⁸O values of authigenic quartz were derived by mass balance, assuming the bulk rock composition can be treated as a three-component mixing problem. We use measured δ¹⁸O of bulk rock, prescribed detrital quartz and clay δ¹⁸O values, typical shale clay/detrital quartz ratios, and measured Al concentrations to estimate the fraction of oxygen in authigenic quartz and its δ¹⁸O value (all data given in SI Appendix, Supporting Information 4). This approach assumes clay minerals—dominated by illite—are the primary Al host, an assumption validated by the XRD spectra (SI Appendix, Supporting Information 1).

We present the most comprehensive chert δ¹⁸O compilation to date (Fig. 1 and SI Appendix, Supporting Information 5). Some anomalously low δ¹⁸O-values in Archean cherts are unlikely to be explained by sediment advection along geothermal gradients but by a different mode of formation. For instance, cherts interbedded in volcanic rocks were subject to locally very high temperatures (2), and narrow microscale δ¹⁸O distributions likely reflect the reequilibration of δ¹⁸O during hydrothermal alteration at 200 to 300 °C (59). In these cases, low δ¹⁸O values do not represent the integrated t–T path during burial diagenesis. We have screened the database to include only sedimentary chert and omit silicified protoliths and vein mineralizations of SiO₂.

We have developed a 1D reaction–advection–diffusion model for silica diagenesis in marine sediments. The model tracks 12 independent mass balances (¹⁶O, ¹⁸O, and Si in four phases: opal-A, opal-CT, quartz, and porewater). Temperatures at each depth are defined from the heat flow and the geometric mean of water and sediment thermal conductivities and are subsequently used to define water–silica oxygen isotope fractionation factors, silica polymorph solubilities, diffusion coefficients, and silica polymorph dissolution and precipitation rate constants. Temperature profiles are dependent on sediment thermal conductivities, porosities, and compaction length scale, so we investigate three depositional settings (margin, shelf, and abyss) using typical values for these parameters in each setting from ref. (60). The model is run for each ensemble of input parameters until a depth is reached at which >99.99% of silica is present as quartz, and δ¹⁸O_quartz has achieved a stable value. Once other local processes affecting porewater δ¹⁸O are accounted for, the model successfully produces the depth of precipitation and δ¹⁸O of authigenic quartz in Sea of Japan ODP cores 795 and 799 (27). Our base model parameterization considers sediment thermal conductivity λ_s = 2.5 W m⁻¹ K⁻¹, initial porosity φ₀ = 0.74, sediment compaction length scale c = 1.7 × 10⁻³ m⁻¹, sedimentation rate  = 100 m Myr⁻¹, bottom water temperature = 25 °C, bottom water dissolved Si concentration = 100 μM L⁻¹, and detrital aluminosilicate mass fraction f_sil = 0.05. Unless specified, all results assume a modern average heat flow of 80 mW m⁻². All further parameters and constants are given in SI Appendix, Supporting Material 2.

M.T., P.J.F., and M.O. designed research; M.T., P.J.F., and M.O. performed research; D.H., N.K.L., and M.W. analyzed data; M.T. collected samples; P.J.F. created model; and M.T. and P.J.F. wrote the paper.

The authors declare no competing interest.