A High‐Resolution Belemnite Geochemical Analysis of Early Cretaceous (Valanginian‐Hauterivian) Environmental and Climatic Perturbations

The Early Cretaceous Weissert event, characterized by a positive carbon isotope excursion and coincident with the Paraná‐Etendeka volcanism, saw a biogeochemical chain of events that ultimately led to an increase in carbon burial. A conclusive link between the Paraná‐Etendeka volcanism and its impact upon the environment remains, however, elusive. Here we reconstruct temperature through the Weissert event from Mg/Ca ratios of belemnites from the Vocontian Trough (France) and SE Spain and use carbon isotopes to link our temperature reconstruction to marine records of carbon cycling. We provide evidence that the Paraná‐Etendeka volcanism, unlike some large igneous provinces, did not cause a climate warming. The case can be made for cooling in the last stages of the Weissert event, which possibly reflects substantial CO2 drawdown. In the absence of warming and consequent accelerated hydrological cycling and the relatively long duration of the eruptive phase of the Paraná‐Etendeka, an alternate trigger for increased fertilization of the oceans is implicated.


Introduction
The Early Cretaceous Weissert event, characterized by a positive carbon isotope excursion, coincides with an increase in atmospheric CO 2 , a demise of shallow-water carbonate platforms (Föllmi et al., 2006), and biocalcification crisis (Barbarin et al., 2012). The origin of the Weissert event has been the subject of much debate whereby Lini et al. (1992), Gröcke et al. (2005), Erba et al. (2004), and Duchamp-Alphonse et al. (2007) have speculated that an increase in atmospheric CO 2 and environmental changes are linked to volcanism associated with the Paraná-Etendeka igneous province. Although links between other large igneous provinces (LIPs) and warming events have been established (e.g., Kerr, 1998;Tejada et al., 2009), this is not the case for the Paraná-Etendeka igneous province (e.g., Bodin et al., 2015;Dodd et al., 2015). Following the positive carbon isotope excursion of the Weissert event, subsequent organic-carbon storage has been considered to have triggered a decline in atmospheric pCO 2 . Many studies (e.g., McArthur et al., 2004McArthur et al., , 2007Meissner et al., 2015;Price & Passey, 2013;Pucéat et al., 2003;van de Schootbrugge et al., 2000) do indeed propose cooling in the late Valanginian (coincident with the most positive carbon isotopes values), based on oxygen-isotope records derived from belemnites and fish tooth enamels. In response to a period of global cooling, McArthur et al. (2007) propose also the formation of a substantial amount of polar ice. The occurrence of glendonites and dropstones (Kemper, 1987;Price, 1999) and changes in calcareous nannofossils (Melinte & Mutterlose, 2001) also appear to confirm cold conditions during the Valanginian. Despite these numerous studies, a robust picture of Valanginian temperature evolution has yet to be fully developed. Here we evaluate high-resolution temporal trends through the Early Cretaceous (Valanginian-Hauterivian) of oxygen isotopes and Mg/Ca ratios of belemnites from the Vocontian Trough (France) and SE Spain in order to develop a robust view of the evolution of temperature. The analysis of carbon isotopes permits the comparison of paleotemperatures with changes in carbon cycling.
During the Valanginian-Hauterivian interval the Subbetic domain was situated on the southern margin of the Iberian plate in the western region of the Tethys Ocean, at a low latitude between 20°and 25°N (Masse et al., 1993;Figure 2). From the Latest Jurassic to the Early Cretaceous, it was a passive margin, with thick hemipelagic postrift sedimentation smothering the submarine relief (Martín-Algarra et al., 1992).

Materials and Methods
The sections studied in Spain, located near Caravaca and Cehegín, are as follows: Prado Borda section, 2.1 km SSE of Caravaca (38°05 0 13°N 1°51 0 14°W); sections on the northeastern slope of Mai Valera Mountain, 2 km west of Cehegín (38°05 0 39°N 1°49 0 09°W); sections along the Cañada Luenga ravine, 3.1 km SSW of Cehegín (38°04 0 02°N 1°48 0 49°W); and the Barranco del Garranchal section, 5.5 km SSW of Cehegín (38°02 0 42°N 1°49 0 20°W). The sections studied in France are the La Charce section (44°28 0 09°N 5°26 0 37°W), an expanded, well-documented section that is the currently suggested Global Boundary Stratotype Section and Point (GSSP) candidate for the base of the Hauterivian Stage (Gradstein et al., 2012); the Vergol/Morénas section (a Berriasian/Valanginian GSSP candidate; Gradstein et al., 2012, 44°12 0 12°N 5°25 0 09°W); and the Angles section (the Valanginian Hypostratotype, 43°56 0 28°N 6°32 0 32°W). Belemnite samples from each section were collected bed by bed, and whenever possible, multiple samples were collected from each bed. In total, 136 belemnites from France have been analyzed, of which 35 of these are excluded due to potential diagenesis (see below) and are combined with a further 226 published analyses (from Bodin et al., 2015;Li et al., 2013;McArthur et al., 2007;van de Schootbrugge et al., 2000). In total, 124 belemnites from Spain have been analyzed in this study, of which 9 of these are excluded due to potential diagenesis and are combined with a further 18 published analyses (from McArthur et al., 2007). In order to accommodate these published isotopic data and isotopic data derived from different sections, data are plotted as numerical ages based on the ages of Tethyan ammonite zones/subzones from Reboulet et al. (2014) in Martinez et al. (2015). Ages are determined according to position of each sample within each ammonite zone/subzone, assuming a constant sedimentation rate. Adjustments for a variable sedimentation rate are made with respect to different thickness of each ammonite zone/subzone at each section examined.
The preservation of the belemnite rostra was assessed using cathodoluminscensce ( Figure 3) using an MK5 CITL instrument and trace element analysis (Ca, Sr, Mg, Fe, and Mn concentrations). The belemnites were prepared for stable isotope and trace element analysis by first removing the areas of the rostrum typically most prone to diagenesis (the rostrum exterior, apical region, alveolus, and observable cracks/fractures; e.g., McArthur et al., 2007;Meissner et al., 2015). The remaining calcite was then fragmented, washed in ultrapure (Type 1) water, and dried in a clean environment. Using 300 to 400 μg of carbonate, stable isotope data were generated on a VG Optima mass spectrometer with a Gilson autosampler at the University of Plymouth. Isotope ratios were calibrated using standards of National Bureau of Standards and were given in δ notation relative to the Vienna Pee Dee Belemnite (V-PDB). Reproducibility was generally better than 0.1‰ for samples and standard materials. The subsamples taken for trace element analysis were digested in HNO 3 and analyzed by inductively coupled plasma-atomic emission spectrometer using a PerkinElmer 3100 at the University of Plymouth. Based upon analysis of duplicate samples reproducibility was better than ±3% of the measured concentration of each element. Repeat analyses of standards Japanese limestone standard-1 and British Chemical Standard Certified Reference Material 393 were within 2% of the certified values for Sr, Mn, Ca, and Mg and 10% for Fe.

Results
Belemnites analyzed from Spain include Duvalia tornajoensis, D. cf. lata constricta, D. binervia, D. cf. emericii, Hibolithes, H. cf. jaculoides, Berriasibelus, Castellanibelus, and Pseudobelus, consistent with earlier studies (e.g., as illustrated in Janssen, 1997Janssen, , 2003. From France, Valanginian and Hauterivian belemnites include Duvalia, Hibolithes, Berriasibelus, Castellanibelus, and Pseudobelus. The belemnite oxygen isotopes ( Figure 4 and supporting information data) show only modest variability ranging from À0.5‰ to 0.3‰ (V-PDB) from Spain and a wider range (À4.0‰ to 2.2‰, V-PDB) from France. The carbon isotope values of the belemnites range from À2.6‰ to 1.9‰ from Spain and À4.9‰ to 1.8‰ from France. Most belemnites sampled in this study were translucent, although some samples from France exhibited prevalent cloudy and opaque areas particularly around the margins of the rostra and the originally partially porous apical region (Figure 3). These areas tended to be Mn rich as revealed bright orange luminescence. As noted above, such areas were either removed prior to or avoided during subsampling. The determined elemental ranges of belemnite rostra (Table S1 in the supporting information) were as follows: Sr (366-2,898 ppm); Mn (1-851 ppm); Mg (1,624-4,684 ppm); and Fe (21-9,186 ppm). The Ca concentrations ranged from 28.7% to 47.8%. Elemental ranges for Spain and France were very consistent. Further, no tendency for a particular belemnite species to show elevated levels of Fe and/or Mn and hence possibly be more susceptible to diagenetic alteration was noted (Figure 3). In line with other data sets from Spain and France (e.g., Armendáriz et al., 2013;Bodin et al., 2009;van de Schootbrugge et al., 2000) those samples with high Fe and Mn concentrations were considered likely to have undergone some isotopic exchange registered by the precipitation of postdepositional diagenetic calcite and were hence excluded from any further analysis. These high Fe and Mn samples are also typically characterized by low Sr concentrations. Indeed, belemnite samples with low Sr concentrations (<800 ppm) may also suggest some diagenetic alteration (van de Schootbrugge et al., 2000;Wierzbowski et al., 2013) and so were also excluded from any further analysis.
The oxygen isotope data from France fluctuate around 0.0‰ ( Figure 5). The most positive values are recorded in the Neocominensiformis Zone, with lighter values seen in the Inostranzewi-Verrucosum and Loryi-Nodosoplicatus Zones. Thereafter, δ 18 O belemnite values fluctuate around À0.4‰. The oxygen isotope data from Spain show a similar pattern with values fluctuating around À0.5‰ (Figure 6). The most positive values are recorded in the Peregrinus Zone. In terms of carbon, a major positive isotope excursion is evident in the data from France and Spain, beginning in the Inostranzewi ammonite zone and culminating in the Verrucosum Zone (Figures 5 and 6). This excursion reveals a 3‰ change in values and is correlated to the Weissert event (e.g., Erba et al., 2004). Bulk rock data (e.g., Duchamp-Alphonse et al., 2007;Hennig et al., 1999;Kujau et al., 2012;Martinez et al., 2015) likewise show the excursion beginning in the Inostranzewi Zone. Also shown in Figures 5 and 6 are Mg/Ca (mmol/mol) ratio derived paleotemperatures. A difference in mean Mg/Ca between Duvalia and Hibolithes of 20% is observed. To test if this was a statistical difference, a Student t test was conducted. The difference was considered significant at P < 0.05. Many studies of seawater temperatures and calcitic Mg/Ca ratios (e.g., Klein et al., 1996;Surge & Lohmann, 2008;Vander Putten et al., 2000) have also revealed interspecific differences in temperature sensitivity or related to taxa living at different depths. We therefore minimize this species-specific effect, by using a normalizing approach (e.g., McArthur et al., 2007) whereby we reduced Mg concentrations in Hibolithes by 20% (Table S1 reports unadjusted data). Mg/Ca paleotemperatures were calculated here using the equation of Lear et al. (2002) for low Mg calcite (benthic foraminifera) and the calcification constants suggested by Bailey et al. (2003) for belemnites: (where Mg/Ca is in mmol/mol) which is consistent with other belemnite studies (e.g., McArthur et al., 2007). However, since the biomineralization of Mg/Ca in belemnite calcite cannot be observed or measured, then the paleotemperature values should be treated with caution. Other Mg/Ca temperature equations (e.g., Klein et al., 1996;Surge & Lohmann, 2008) provide different absolute temperatures. For example, applying the Mg/Ca-temperature calibration for the extant oyster Crassostrea virginica (Surge & Lohmann, 2008) to our fossil belemnite data results in temperatures that are considerably lower (by~4-8°C) than those derived using the equation of Lear et al. (2002). Cretaceous seawater Mg/Ca was likely to be significantly lower compared to modern seawater (e.g., Stanley & Hardie, 1998), and this could also be a significant source of uncertainty and may have led to Cretaceous calcifiers being characterized by lower Mg/Ca ratios. Hence, only Mg/Ca-derived paleotemperature trends (and derived δ 18 O seawater trends; see below) can be examined with confidence. Figure 6. Record of oxygen isotopes, carbon isotopes, Mg/Ca temperatures, and δ 18 O seawater (with Locally weighted smoothing using the PAST software package (Hammer et al., 2001), through the Valanginian interval from Spain plotted against numerical age for the Valanginian. Blue curves indicate the 95% confidence interval derived using a bootstrap technique. Numerical calibration is modeled on sediment thickness, with adjustments for a variable sedimentation rate. Numerical ages based on using boundary ages of Martinez et al. (2015). Tethyan ammonite zones from Reboulet et al. (2014). The plots incorporate the belemnite isotope and Mg/Ca data of McArthur et al. (2007), plotted as open circles. V-PDB = Vienna Pee Dee Belemnite.
In both France and Spain temperatures using this methodology, modest variability is observed during the Berriasian and earliest Valanginian, and a pronounced cooling event is recorded across the Verrucosum/Peregrinus zonal boundary following the most positive carbon isotope values of the Weissert event. Using the calculated paleotemperatures from Mg/Ca, we then use the equation of Kim and O'Neil (1997) to derive δ 18 O seawater (Figures 5 and 6). Utilizing this approach, mean values of the oxygen isotopic composition of seawater from France are calculated to fluctuate around 1.0‰ with more positive values seen within the Neocominensiformis and Furcillata-Loryi Zones. More negative δ 18 O seawater values are notable in the Verrucosum Zone. Lighter values are also calculated for the Verrucosum Zone data from Spain.

Discussion
A major positive carbon isotope excursion is evident in the data from France and Spain, documenting a significant perturbation of the carbon cycle (Figures 5 and 6). A relatively high degree of scatter is also observed. Different belemnite species may have recorded different conditions in the water column, related to differences ecology. However, Duvalia and Hibolithes (the two most commonly encountered species) do show overlapping ranges of carbon isotope values ( Figure 4 and Table S1). This suggests a limited species/habitat effect upon carbon isotopes during biomineralization. Alternatively, it might be assumed that some altered samples may not have been identified by the trace element and petrological screening and isotopic outliers result from this. Nevertheless, this carbon isotope excursion appears correlatable to the Weissert event.
The volcanic activity associated with emplacement of the continental Paraná-Etendeka LIP has been considered as a trigger of environmental change associated with the Valanginian Weissert event (e.g., Bajnai et al., 2017;Charbonnier et al., 2017;Duchamp-Alphonse et al., 2007;Erba et al., 2004;Lini et al., 1992;Martinez et al., 2015;Weissert et al., 1998). It might be expected that δ 13 C-depleted volcanic CO 2 (with a δ 13 C of approximately À5‰; Kump & Arthur, 1999) from the Paraná-Etendeka would have resulted in the ocean/atmosphere recording a negative event. This is clearly not the case (Figures 5 and 6), as the Weissert event is characterized by a positive carbon isotope excursion. Notably, the much larger Siberian Trap volcanism, which straddles the Permian-Triassic boundary (e.g., Reichow et al., 2002), is associated with a prominent negative carbon-isotope excursion, but volcanic CO 2 is also unlikely to have been a main cause of the negative excursion (Kump & Arthur, 1999). A reevaluation of the timing of the onset of the Weissert event is ascertained to be 135.22 ± 1 Ma derived from U-Pb ages from tuff layers in the Neuquén Basin (Aguirre-Urreta et al., 2015) and an update of the Valanginian-Hauterivian astrochronological time scale (Martinez et al., 2015). Therefore, the Valanginian Weissert event appears to have coincided with the onset of the eruptive phase of the Paraná-Etendeka, which has recently been dated between 134.6 ± 0.6 Ma and 134.3 ± 0.8 Ma (Janasi et al., 2011;Thiede & Vasconcelos, 2010). Charbonnier et al. (2017) also show mercury enrichments linked to Valanginian volcanism within the Pertransiens-Inostranzewi Zones. Dodd et al. (2015) recently presented a detailed magnetostratigraphy for the Etendeka portion of the province that suggested emplacement took place a little earlier during Chron 15 (Figures 5 and 6) and suggest a minimum period of volcanic activity in excess of 4 Myr.
The volcanism of the Paraná-Etendeka presumably increased atmospheric CO 2 , which in turn should elicit a warming. The oxygen isotope values from both Spain and France show very little response to the coincident Paraná-Etendeka LIP (Figures 5 and 6). Either the anticipated warming is erroneous or that the oxygen isotope values are not recording the warming and other factors such as the isotopic composition of seawater is compromising the projected oxygen isotope-derived temperature record. Mg/Ca ratios provide a methodology that is free from the effects of changes in the oxygen isotopic variability of seawater but sensitive to other paleoenvironmental drivers. The Mg/Ca ratio in seawater is spatially constant and unlikely to change on timescales of less than 1 million years due to the long residence times of both Mg and Ca in the oceans (Broecker & Peng, 1982), although, as noted above, Cretaceous seawater Mg/Ca was likely to be lower than modern seawater (e.g., Stanley & Hardie, 1998). Laboratory culture studies have shown that the pH, carbonate ion concentration, and salinity of seawater act as controls on Mg/Ca ratios in foraminifera but suggest that their influence is small in comparison with temperature (Elderfield et al., 2006;Lea et al., 1999). More recently some studies (e.g., Dueñas-Bohórquez et al., 2009;Ferguson et al., 2008) have suggested that salinity may have had a greater than hitherto expected influence on foraminiferal Mg/Ca ratios. Dueñas-Bohórquez et al.
(2009) suggest that an increase of four salinity units is equivalent to about 1°C warming, in terms of Mg/Ca ratios. Hence, the changes observed in the Mg/Ca data (Figures 5 and 6) are considered to be temperature related, although a minor contribution from salinity/evaporation cannot be excluded.
The Mg/Ca-derived temperatures from Spain show little variability during the Berriasian and earliest Valanginian. From France, a similar pattern is observed with a cooling event seen during the Neocominensiformis ammonite zone. A coeval cooling event is not present in Spain and could relate to the lower-resolution data of the Spain section or be a localized event. In both records (Figures 5 and 6) following an interval of warmth, a pronounced cooling event is recorded across the Verrucosum/Peregrinus zonal boundary following the most positive carbon isotope values of the Weissert event. The longer record from France sees a recovery of temperatures during the latest Valanginian and Hauterivian. Hence, unlike other intervals marked by LIP eruptions (e.g., the Ontong Java Plateau, Tejada et al., 2009, or the Caribbean, Kerr, 1998, it is clear from the data presented here that there are no substantial temperature changes (warming) occurring during the main projected phase of the Paraná-Etendeka LIP, despite the large area and volume (approximately equivalent to the Siberian traps, 1-2 × 10 6 km 3 ; Wignall, 2001). It is conjectured that the longer duration of volcanism (4-5 Myr), compared to other LIPs of a similar volume, increased atmospheric and biological recovery time between individual eruptions (Dodd et al., 2015).
A case can be made for cooling in the last stages of the Weissert event ( Figures 5 and 6), which possibly reflects substantial CO 2 drawdown. Indeed, increased fertilization of the oceans (e.g., Duchamp-Alphonse et al., 2007;Erba et al., 2004) causing sequestration of marine organic carbon and a consequence decrease (drawdown) of atmospheric CO 2 has been suggested. In the absence of warming and consequent accelerated hydrological cycling (cf. Erba et al., 2004), and the relatively long duration of the eruptive phase of the Paraná-Etendeka, an alternate trigger for increased fertilization of the oceans is implicated. For example, a coupling of cooling, higher dust flux, and increased iron fertilization and productivity has been observed, for example, during the Quaternary (e.g., Martínez-García et al., 2014).
A temperature recovery occurs after the Paraná-Etendeka episode and the Weissert event, to pre-Weissert temperature levels in the latest Valanginian. The inferred cooling coincides with a stable or decrease in the estimated δ 18 O seawater . Hence, it is difficult to envisage substantial amounts of ice, as in such circumstances a trend to more positive δ 18 O seawater values (globally) would be expected (cf. McArthur et al., 2007). The Valanginian-Hauterivian boundary does see positive δ 18 O seawater values, and following the above reasoning, δ 18 O seawater values could imply polar ice locking up light oxygen isotopes. However, a cooling trend in the Mg/Ca temperature data is not apparent; hence, δ 18 O seawater trend could be related to regional changes , such as changes in evaporation or ocean circulation.

Conclusions
This study evaluates the links between the Paraná-Etendeka volcanism and the Weissert event. Unlike some LIPs, the data presented here indicate that the Paraná-Etendeka did not cause a climate warming (or a mass extinction; Wignall, 2001). The case can be made for cooling in the last stages of the Weissert event, which possibly reflects substantial CO 2 drawdown. In the absence of warming and consequent accelerated hydrological cycling (cf. Erba et al., 2004), and the relatively long duration of the eruptive phase of the Paraná-Etendeka, an alternate trigger for increased fertilization of the oceans is implicated. We interpret the δ 18 O belemnite signal as being primarily driven by δ 18 O seawater changes and not paleotemperature.