Stepwise oxygenation of the Paleozoic atmosphere

Oxygen is essential for animal life, and while geochemical proxies have been instrumental in determining the broad evolutionary history of oxygen on Earth, much of our insight into Phanerozoic oxygen comes from biogeochemical modelling. The GEOCARBSULF model utilizes carbon and sulphur isotope records to produce the most detailed history of Phanerozoic atmospheric O2 currently available. However, its predictions for the Paleozoic disagree with geochemical proxies, and with non-isotope modelling. Here we show that GEOCARBSULF oversimplifies the geochemistry of sulphur isotope fractionation, returning unrealistic values for the O2 sourced from pyrite burial when oxygen is low. We rebuild the model from first principles, utilizing an improved numerical scheme, the latest carbon isotope data, and we replace the sulphur cycle equations in line with forwards modelling approaches. Our new model, GEOCARBSULFOR, produces a revised, highly-detailed prediction for Phanerozoic O2 that is consistent with available proxy data, and independently supports a Paleozoic Oxygenation Event, which likely contributed to the observed radiation of complex, diverse fauna at this time.

for consistency, however, using the original GEOCARBSULF value of 35.2‰ and the highest fractionation (48.1‰) produces a δ 34 S record which for the last 250 Ma is near identical to the geologic record.

Supplementary Figure 1 | Ocean sulphate isotope predictions for the Phanerozoic from
our GEOCARBSULFOR model. The red line is the ocean δ 34 S output using a static average isotopic fractionation from 570 Ma to Present (35.95‰), calculated from the Wu et al. 1 dataset. The grey envelope is the ocean δ 34 S output using the lowest (16.7‰) and highest (48.1‰) fractionation values from Wu et al. 1 The black line is the geologic δ 34 S data from Wu et al. 1 We conducted a test to ascertain whether simply substituting the sulphur isotope fractionation equation in the original GEOCARBSULF with the Wu et al. 1 fractionation data, would produce low pO2 in the early Paleozoic, as Wu et al. 1 take into account the effects of reoxidation and disproportionation on sulphide (though not the possibility of other sinks, such as organic sulphur) changing the isotopic signature of sedimentary pyrite. This change leads to unfeasibly high (considering geochemical redox data e.g.2,3 ) oxygen concentrations of >35% atm in the early Paleozoic, and negative O2 during the late-Mesozoic to Cenozoic. Therefore, a pO2 dependent feedback embedded somewhere in the sulphur cycle is a requirement for the model to produce reasonable levels of atmospheric O2.
As GEOCARBSULFOR can calculate the amount of pyrite and sulphate being buried at each time-step, we examined whether the model was generating reasonable values for the fraction of total sulphur leaving the ocean as pyrite (fpyr), by comparing the model outputs to other modelling approaches 4,5 . We tried two different methods to compute fpyr: one was a direct method, using the burial fluxes, following supplementary equation (57) 5 . As with our synthetic ocean δ 34 S sulphate record in Supplementary Fig. 1, these fpyr results suggest a transition from lower to higher (persistently >35.9‰the average value over the last 570 Ma) sulphur isotope fractionation values in the Permian, and this can be seen in the geologic record. In GEOCARBSULF 6 , Berner contends that erosive stripping is a more important process in the weathering of pyrite (and organic matter) than oxygen levels, and although the recent iteration of COPSE 7 concurs, the original COPSE 8 model and others 9,10 , include a dependency on atmospheric oxygen levels when considering the weathering of reductants (e.g. pyrite). In Supplementary Fig. 3 we conduct a sensitivity analysis of the model to an oxidative feedback term in the equations for young and ancient pyrite weathering, as well as an analysis of the importance of linking young pyrite weathering to the weathering of silicates.
In Supplementary Fig. 3, the black line is our GEOCARBSULFOR baseline, while the grey line is almost the same, but with the oxidative feedback (O2mr 0.5 ) removed from the pyrite weathering equations (equation (3) and supplementary equation (17)). The magenta line is a model identical to GEOCARBSULFOR but the pyrite weathering equations are those used in the original GEOCARBSULF 6 , while the green line is the same as the magenta line, but with oxidative feedback included into Berner's equations. Finally, the yellow line is as the magenta line, but with a fixed uplift rate of 1, for the entire run.
Our results show that although there is some sensitivity to the inclusion of an oxidative feedback term in pyrite weathering, it does not alter the general trend of the evolution of pO2 as predicted by GEOCARBSULFOR. Comparison of our baseline GEOCARBSULFOR Two final tests were conducted, to assess the robustness of our model outputs. These tests were to evaluate both the amount of sulphate assigned to the ocean-atmosphere reservoir by the model (Supplementary Fig. 4), and the CO2 levels generated ( Supplementary Fig. 5), against other models and geochemical proxies. In both tests, GEOCARBSULFOR outputs compare favourably with both the proxy data and other modelling efforts. CO2 levels predicted by GEOCARBSULFOR are on the low side of the proxy bands, particularly during the Cretaceous and Paleogene. The Paleogene low is a common issue with previous work involving GEOCARBSULF [11][12][13] , and requires the seafloor spreading rate to be ~3x the present day value in order to reconcile the model predictions with the proxy record. The discrepancy between model predictions and proxy estimates is possibly due to underrepresentation, or absence, of all plate boundary types and their contribution to CO2 degassing 14 . However, recent work on this issue has resulted in GEOCARBSULF predictions that better match the Paleogene record, but are subsequently too high for the Neogene 14 .
Nevertheless, we believe that our model has been comprehensively checked against the available geologic data and provides a good fit to our understanding of various geochemical indices, and thus GEOCARBSULFOR gives a compelling indication of low pO2 in the early Paleozoic.

Supplementary Note 2: Full model parameters
For the following tables, where we mention GEOCARBSULF as a source, we mean the final values as used in the 2009 22 paper, some of which are unchanged from, for example, GEOCARB II 23 .  For all time arrays we apply interpolation between data points.

Supplementary Note 3: Full model equations
The full set of equations for the model are documented below, except the new burial equations for pyrite and gypsum, and weathering equations for young pyrite and young gypsum (equations 3 -6), which can be found in the Methods section of the main paper.
These are taken from the COPSE 8 model. All other equations have been taken from GEOCARBSULF, but we note where alterations to the original equations have been made.