Age model and study sites

Samples from DSDP Site 588 and ODP sites 761 Hole B, 926 Holes A/B, and 872 Hole C were used to reconstruct SST and \(\left[ {{\mathrm{B(OH)}}_4^ - /{\mathrm{DIC}}} \right]_{{\mathrm{sw}}}\) across the early to middle Miocene (22.0–11.5 Ma; Supplementary Fig. 1). Planktic foraminifera Mg/Ca, B/Ca, and δ13C p records from ODP site 761B were sampled at an average temporal resolution of 23 kyr. In contrast, our other sites were sampled at lower resolution for B/Ca: ODP 926 (~200 kyr) and DSDP 588 (~100 kyr) and ODP 872 (~1 Myr). ODP 761 (16°44.23’S, 115°32.10’E, water depth of 2179 m) is situated in the Indian Ocean on the Wombat Plateau. Site ODP 926 (3°43.148’N, 42°54.507’E, water depth 3598 m) in the eastern equatorial Atlantic Ocean on Ceara Rise. Site DSDP 588 A (26°06.7’ S, 161°13.6’ E, water depth 1548 m) in the Southwest Pacific Ocean on Lord Howe Rise. ODP Site 872 (10°N, 162°E, 1287 m water depth) is in the tropical north Pacific gyre on the sedimentary caps of flat-topped seamounts.

The age model for Site 761B is based on a fourth-order polynomial fit through the biostratigraphic and isotopic datums14,35. For ODP site 926, we use a polynomial fit through nannofossil and planktic foraminifer biostratigraphic datums adapted to the CK95 timescale35. The published age model for DSDP 58859 is based on compiled biostratigraphy, magnetostratigraphy, isotope stratigraphy calibrated to CK95. Ages for ODP Site 872 were calculated by linear interpolation between reliable biostratigraphic datums60 from the early to mid-Miocene adapted to the CK95 timescale35. Core sites29,30,35,61 used in this study were selected because currently they are located in regions where modern surface waters are close to equilibrium with respect to atmospheric CO 2 (Supplementary Fig. 1). We use δ11B-pH estimates35 from ODP sites 761, 872, and 926 in addition to B/Ca presented in this study to estimate surface ocean \(\left[ {{\mathrm{B(OH)}}_4^ - /{\mathrm{DIC}}} \right]_{{\mathrm{sw}}}\) and\(\left[ {{\mathrm{B(OH)}}_4^ - /{\mathrm{HCO}}_3 - } \right]_{{\mathrm{sw}}}\).

Planktic foraminiferal taxonomy and ecology

We generated trace element and isotope records using T. trilobus, a multi-chambered, photosymbiont bearing species, which is predominantly a mixed layer dweller calcifying at 0–50 m and is abundant in subtropical and tropical oceans. T. trilobus is a morphospecies of Trilobatus sacculifer62,63 present throughout the Neogene and Quaternary periods. T. trilobus has been used extensively in previous Miocene climate and boron isotope studies29,23,30,35 due to its narrow habitat range, well-defined δ11B-pH calibration, calcification close to equilibrium conditions (i.e., minimal vital effects), and abundance29,35,63,64. Further, previous studies demonstrated that temperatures derived from T. sacculifer are most suitable for estimating annual mean SST65 in tropical waters (between 20° N/S) within ±1 °C. Visual inspection of preservation (Supplementary Fig. 4) of T. trilobus from ODP site 761B show moderately good preservation with no visible signs of infilling or dissolution. Visual inspection of foraminiferal specimens from ODP site 926 shows evidence of a recrystallized nature, however the Mg/Ca and boron isotope signature has remained intact29,35,66. Visual inspection of foraminiferal specimens from ODP site 872 show a glassy nature indicative of well-preserved specimens and exceptionally good carbonate preservation supported by previous work using scanning electron microscopic (SEM) images60,67. Visual inspection of foraminiferal specimens from DSDP 588 have been shown to be pristine in nature with no evidence of dissolution and have been previously used in oxygen isotope reconstructions52.

Trace metal and isotopic analysis

Between 30–40 tests of the planktic foraminifera T. trilobus were picked from the 300–355 μm size fraction from ODP site 761, 926, and DSDP 588. ODP 872 samples were prepared for boron isotope analysis and an aliquot of the foraminiferal sample was used to collect trace element data (see ref. 35). Picked specimens were weighed, gently crushed between glass plates and homogenized for chemical cleaning and geochemical analysis. In samples where T. trilobus abundance was low (~30 samples), fewer specimens (10–20 individuals) were analyzed. After sample homogenization, an aliquot of the foraminiferal sample was used to collect trace element and δ13C p data. It is noteworthy that only δ13C p from ODP site 761 are presented here.

Test fragments for Mg/Ca and B/Ca analyses were cleaned using a protocol to remove surficial clays, adhered organic matter, and potential secondary carbonate overgrowths68. Clay removal consists of repeated rinses with low boron Milli-Q water and methanol. The oxidative cleaning step consists of a sodium hydroxide and hydrogen peroxide solution to remove any present organic matter. The reductive cleaning step was not included in the cleaning process as this step is unnecessary for B/Ca as the cleaning does not alter the B/Ca ratio of foraminifera or calcite bound boron69. Lastly, a weak acid leach using 0.001 M hydrochloric acid was conducted to remove any ions re-absorbed during the cleaning process on the test surface. Following the clay removal and oxidative steps, samples were examined under a binocular microscope and visible non-carbonate particles were removed using a fine paintbrush. Cleaned, treated samples were dissolved in trace metal pure 0.065 M HNO 3 and diluted with trace metal pure 0.5 M HNO 3 to a final volume of 350 μl, to achieve a target calcium concentration of 4 mM. Trace element analysis (B, Mg, Al, Ca, Mn, Fe) for ODP 761, ODP 926, and DSDP 588 was carried out at Cardiff University on a Thermo Scientific Element XR Sector Field Inductively Coupled Plasma Mass Spectrometer (SF-ICP-MS), whereas samples from ODP site 872 were analyzed at the University of Southampton, using a Thermo Scientific Element 2 XR SF-ICP-MS. In both instances, calcium concentrations of bracketing standards were matched to foraminifera samples to reduce matrix effects70,71.

To monitor the possibility of clay contamination Mg/Ca and B/Ca data were rejected if Al/Ca exceeded 80 μmol/mol for all sites. Additional cleaning effectiveness was supported by observing no significant correlation between Mg/Ca, and Fe/Ca and Mn/Ca (Supplementary Figs. 5-6). There is no correlation between B/Ca and Fe/Ca at any of the investigated sites. There is however a positive correlation between Mn/Ca and B/Ca at ODP site 926 (R2 = 0.4) and 761 (R2 = 0.6), and a slight negative correlation (R2 = 0.3) between Mn/Ca and B/Ca at ODP site 872 (Supplementary Fig. 6). However, at site 761, correlations across 500 kyr windows from 16.5 to 11.5 Ma show no significant correlation (R2 < 0.3; p > 0.05) between these variables in 8 out of the 10 intervals (Supplementary Table 1). We therefore propose here that the Mn/Ca record reflects another aspect of the system (e.g., oxygenation, productivity), which would be correlated on geological timescales with B/Ca as in the Miocene72. Long-term precision at Cardiff University and University of Southampton was determined by analyzing an independent consistency standard during each run for 1 year. Reported values are 0.5% and 4% (r.s.d.) for Mg/Ca and B/Ca, respectively, for analyses at Cardiff University and 2% and 4% (r.s.d.) for Mg/Ca and B/Ca, respectively, for analyses at University of Southampton.

Stable carbon isotope ratios were measured at Cardiff University on a Finnigan MAT 252 micro-mass spectrometer Kiel III Carbonate Device when sample weights were <100 μg and measured on a Delta isotope ratio mass spectrometer when samples were greater than 100 ug. Long-term precision based on replicate measurements of a laboratory standard (NBS 19) are 0.08‰ for δ13C.

Mg/Ca-paleotemperature calculation and constraining Mg/Ca seawater changes

Planktic T. trilobus is the morphotype of the modern T. sacculifer, we consider it appropriate to calculate mean annual SSTs from T. trilobus Mg/Ca data using the T. sacculifer calibration65:

$$\frac{{Mg}}{{Ca}} = 0.347\ \left( {\!\pm \! 0.011} \right)e^{0.090\left( {\! \pm \! 0.003} \right)T}$$ (1)

Variations in Miocene foraminiferal Mg/Ca are driven by changes in temperature; however, on longer timescales, both changes in temperature and in seawater Mg/Ca must be considered. Recently, it has been shown that a power function73 best describes this relationship.

$$\frac{{Mg}}{{Ca}}_{foram} = \left[ {\frac{{\frac{{Mg}}{{Ca}}_{sw}\left( t \right)}}{{\frac{{Mg}}{{Ca}}_{sw}\left( 0 \right)}}} \right]^CBe^{AT}$$ (2)

where Mg/Ca sw (t) and Mg/Ca sw (0) are seawater Mg/Ca ratios for the Miocene and present, and A, B, and C are constants (A = exponential, B = pre-exponential, C = power constant), respectively.

Mg and Ca have relatively long residence times (~13 Myr and ~1.1 Myr, respectively) in the ocean. The modern day seawater Mg/Ca value is 5.2 mmol/mol and low-resolution fluid inclusions data74 show a rise across the Neogene to the modern values. Here we use the fluid inclusion value of 3.43 mmol/mol in the paleotemperature calculation. We use the T. sacculifer calibration65 (A = 0.09, B = 0.347) for T. trilobus and apply a power constant75 of 0.41 for T. Sacculifer based on the available data for this species76. T. trilobus SST estimates can be calculated from the following equation:

$$\frac{{Mg}}{{Ca}}_{foram} = 0.293\,e^{0.090T}$$ (3)

Using a range of Mg/Ca seawater estimates, calculated SST during the MCO period in Site 761 ranges from 1.4 °C warmer to uncorrected Mg/Ca-SST estimates being cooler by 0.5 °C than the modern (modern mean annual temperature is ~27.7 °C).

In addition to changes in seawater Mg/Ca, Mg incorporation in foraminifera has shown a species-specific dependency to changes in carbonate chemistry. A recent review77 proposed a pH correction was necessary to correct Mg/Ca records when estimating SST for some species. This study showed that Mg/Ca ratios in T. sacculifer are insensitive to changes in pH. As Mg/Ca ratios are derived from T. trilobus shells, the morphotype of T. sacculifer, we apply no pH correction in this study.

Foraminiferal preservation

Several lines of evidence suggest that dissolution does not significantly affect the Mg/Ca or B/Ca values at Site 761 and suggest our data represent a climate signal. First, ODP Site 761 is situated well above the modern lysocline, above the critical 20 μmol/kg ΔCO 3 2−, in a relatively shallow burial depth during the middle Miocene (<50 m) (GLODAP78). Average shell weight of T. trilobus, from the 300–355 μm size fraction, does not covary with the Mg/Ca or B/Ca record, supporting our argument that these values are not biased (Supplementary Fig. 7).

Despite the reasonable appearance of foraminifera, all tests appear frosty or opaque in contrast to exceptionally well-preserved translucent test shells from hemipelagic muds79. However, this preservation state is typical of most deep-sea carbonates and is caused by micro-recrystallization of calcite. Large-scale recrystallization is not evident in SEM images for ODP site 761 (Supplementary Fig. 4), suggesting that diagenesis did not drive prominent shifts in δ13C and trace elements. In addition, Sr/Ca ratios from ODP site 761 show high values (~1.1–1.2 mmol/mol; Supplementary Table 2), which are consistent across much of the record, suggesting that diagenesis did not have a major influence on Mg/Ca or B/Ca variations. Nonetheless, Mg/Ca values decrease by a negligible amount with initial diagenetic alteration79, thus temporal changes in Mg/Ca are less likely to be affected. Further work80 has shown a decrease in planktic B/Ca in recrystallized relative to well-preserved foraminifera; however, additional work is needed to isolate the primary diagenetic signal. For the reasons outlined above, we believe that diagenesis had a minimal effect on our records and we therefore interpret the geochemical records in terms of paleoceanographic conditions.

Estimation of surface ocean [B(OH) 4 −/HCO 3 −] sw and [B(OH) 4 -/DIC] sw from B/Ca

Similar to boron isotopes, the boron content (expressed as B/Ca) of planktic foraminifera has been suggested to be a proxy for ocean carbonate chemistry69,81. Boron exists in seawater primarily as two species borate (B(OH) 4 − and boric acid (B(OH) 3 ), and the relative concentration of each boron species and their δ11B composition varies with pH in seawater. Boron isotopic evidence and other arguments suggest that B(OH) 4 − is the species predominantly incorporated into the foraminiferal calcite lattice82. B/Ca in planktic foraminifera is controlled by a combination of ocean carbonate chemistry parameters (pH, B(OH) 4 − \({\mathrm{HCO}}_3^ -\) and DIC) as demonstrated through empirical culture experiments and field studies71,81,83. Culturing efforts aimed to disentangle this co-varying carbonate system parameters show that B/Ca is governed principally by the ratio of seawater B(OH) 4 − to DIC or \({\mathrm{HCO}}_3^ -\)81,84. Application of the B/Ca proxy to reconstruct shifts in the concentration of [B(OH) 4 −/DIC] in seawater across major climate transitions has shed light on past carbon cycle perturbations84,85 (e.g., Paleoceane-Eocene Thermal Maximum), where large shifts in the ocean carbonate system occurred. Similarly, in our study we aim to examine the carbon cycle perturbations of the early to middle Miocene, a time interval whereby a large carbonate system shifts likely occurred, suggesting some of the complexities that others have discussed with regard to the glacial–interglacial cycles may not be significant.

Culture experiments with living T. sacculifer demonstrate a positive relationship between B/Ca and oceanic carbonate system parameters, where B/Ca increases at higher pH and lower DIC and HCO 3 − concentrations (Supplementary Fig. 8)81. These experiments also provide calibrations that relate B/Ca to both [B(OH) 4 −/HCO 3 −] sw and [B(OH) 4 −/DIC] sw ; however, they argue for a better fit to the latter parameter. Recent culture work further explored the mechanistic understanding of B uptake, indicating that planktic B/Ca, Orbulina universa specifically, is driven by HCO 3 − in support of [B(OH) 4 −/HCO 3 −] sw rather than [B(OH) 4 −/DIC] sw controlling B/Ca83. Here we present [B(OH) 4 −/HCO 3 –] sw to evaluate the main drivers of carbonate system change through the MCIE. Further, we consider both B/Ca carbonate chemistry sensitivities (e.g., [B(OH) 4 −/HCO 3 −] sw and [B(OH) 4 −/DIC] sw ) to estimate surface ocean DIC.

When we compile all available B/Ca-[B(OH) 4 −/DIC] sw and [B(OH) 4 −/HCO 3 −] sw T. sacculifer calibration data, an offset is observed between culture81 and core-top86 datasets, where the core-top data sits lower and with a slightly steeper slope predicted by culture experiments, although this slope difference may largely be a function of the different ranges spanned by the datasets, as the available core-top data are limited by the range in surface ocean carbonate chemistry (e.g., [B(OH) 4 −/DIC] sw ) that span the depth habitat of wild-type T. sacculifer compared with conditions manipulated in culture work81. We account for differences in test size87 and normalize values to a common salinity (S = 35)81, but an offset and slope difference between the datasets remain (Supplementary Fig. 8). Here we suggest the offset could be related to the fact core-top foraminiferal specimens, during their life cycle, sink to water with lower pH and higher DIC/ HCO 3 − in contrast to the cultured specimens. This is supported by the lower core-top B/Ca values relative to culture. Furthermore, the core-top study used T. sacculifer specimens with the final sac-life chamber, which have lower B/Ca relative to non-sac chambers81, in contrast to the cultured specimens.

To address the B/Ca offset between datasets we retain the slope predicted in culture and normalize the culture to core-top values to generate new ad hoc B/Ca to [B(OH) 4 −/DIC] sw and [B(OH) 4 −/HCO 3 −] sw calibrations (B/Ca = 33 + 773*[B(OH) 4 −/DIC] sw ;B/Ca = 39 + 561*[B(OH) 4 −/HCO 3 −] sw ; Supplementary Fig. 8). Application of these ad hoc calibrations to ODP 761 B/Ca reconstructions requires additional consideration of past changes in seawater boron ([B sw ]) on calibration sensitivity. We can assume a constant seawater composition over the studied interval as it is considerably shorter than the oceanic residence time of B (10–20 million years)88. A positive relationship between seawater B concentration ([B sw ]) and B content in planktic foraminiferal calcite (O. universa) is observed in culture experiments81,83,84 (Supplementary Fig. 9). Following Mg/Ca convention and the observation that Mg/Ca sensitivity changes as a function of seawater Mg/Ca ratio, we adjust the calibration sensitivity to changes in B/Ca sw , to account for the influence of variable Mg/Ca sw on Mg/Ca thermometry70,75 by assuming a linear scaling:

$$\frac{B}{{C{\mathrm{{a}}}}}{\mathrm{{foram}}} = \frac{{\frac{B}{{C{\mathrm{{a}}}}}sw\,t = t}}{{\frac{B}{{C{\mathrm{{a}}}}}sw\,t = 0}}x\frac{{B({\mathrm{{OH}}})_4^ - }}{{{\mathrm{{HCO}}}_3^ - }} \ast {m} + {b}$$ (4)

$$\frac{B}{{C{\mathrm{{a}}}}}{\mathrm{{foram}}} = \frac{{\frac{B}{{C{\mathrm{{a}}}}}sw\,t = t}}{{\frac{B}{{C{\mathrm{{a}}}}}sw\,t = 0}}x\frac{{B({\mathrm{{OH}}})_4^ - }}{{{\mathrm{{DIC}}}}} \ast {m} + {b}$$ (5)

where B/Ca sw (t = t) for the age of the sample and B/Ca sw (t = 0) is modern B/Ca sw and m and b are the slope and y-intercept for each calibration from. Miocene estimates of B sw and Ca sw are derived from boron isotope modeling and fluid inclusion estimates, respectively74,88.

Calculation of DIC and HCO 3 −from the B/Ca-[B(OH) 4 −/DIC] sw and [B(OH) 4 −/HCO 3 −] sw calibrations requires an independent estimate of B(OH) 4 − across the Miocene, and here we use δ11B-pH estimates from ODP 761, 872, and 92635 to calculate B(OH) 4 −. B(OH) 4 − is related to pH as described here:

$${\mathrm{{pH}}} = - \log [H^ + ]$$ (6)

$$B({\mathrm{{OH}}})_4^ - = B_{{\mathrm{{TOT}}}}/\left( {1 + \frac{{[H^ + ]}}{{K_{\mathrm{{B}}}^ \ast }}} \right)$$ (7)

where \(K_{\mathrm{{B}}}^ \ast\) is the stoichiometric equilibrium constant for boric acid89 at time equivalent T, S, and P, and Mg/Ca sw (t = X Ma). pH estimates used were derived from the three δ11B sw scenarios35 (G17, RH13, L02). The pH derived from each scenario have a similar structure but differ in absolute pH values. This should minimally impact the DIC trend across the Miocene, however, we use pH estimates from all scenarios to estimate the likely DIC (Supplementary Fig. 10). Surface ocean pH estimates from ODP 872 are higher by 0.1 pH units compared with site 926 and 761 in the early Miocene (16.5–17.0 Ma) (Supplementary Fig. 2). To account for this offset, we adjust ODP site 872 by 0.1 pH units and use these adjusted pH values in the DIC estimates. SST estimates are derived from Mg/Ca in complementary samples35 and modern salinity is assumed for the entire duration of the Miocene record. We normalize planktic B/Ca values to a common salinity, assuming modern salinity values for each site (S = 35)81. B sw varied in time according to marine boron isotope budget88 in the same manner as above for B/Ca sw estimates. Using the B/Ca-[B(OH) 4 −/HCO 3 −] sw calibration, we estimate DIC from HCO 3 − using the following equation:

$${\mathrm{{DIC}}} = \frac{{[{\mathrm{{H}}}^ + ][{\mathrm{{HCO}}}_3^ - ]}}{{K_1}} + \left[ {{\mathrm{{HCO}}}_3^ - } \right] + \frac{{K_2[{\mathrm{{HCO}}}_3^ - ]}}{{[{\mathrm{{H}}}^ + ]}}$$ (8)

where K 1 and K 2 are the carbonic acid dissociation constants89.

The principal uncertainties in calculating [DIC] sw and [HCO 3 −] in this way were judged to be uncertainties in the B/Ca calibration and in the pH used to calculate [B(OH) 4 −]. Here we used Monte Carlo method to propagate these uncertainties into our final estimates of DIC and HCO 3 −. In order to highlight the long-term trends in DIC and HCO 3 −, we then fit smoothing splines through the data with the degree of smoothing determined by generalized cross validation. Again, the uncertainty in these smoothed trends was determined using a Monte Carlo approach and a consideration of the error bars of each individual estimate.

DIC estimates derived using the range of pH scenarios (G17 v RH13 v L02) show a similar trend across the middle Miocene (Supplementary Fig. 11). Miocene estimates of surface ocean DIC derived from ODP site 761 using both B/Ca-[B(OH) 4 −/DIC] sw and [B(OH) 4 −/HCO 3 −] sw calibrations, and pH estimates from all scenarios show higher DIC during the MCO (14.0–17.0 Ma; Supplementary Fig. 11 and Supplementary Tables 3 and 4; 68% CI; n = 15). Early Miocene (17.0–22.0 Ma) estimates are lower than the MCO (Supplementary Table 3; 68% CI; n = 7) increasing towards the onset of the MCO. MCO DIC estimates decline following the MMCT. MMCT DIC levels (13.5–11.0 Ma) are lower than the MCO (Supplementary Table 3; 68% Confidence Interval; n = 12). Further, the overall decline in DIC from the MCO to post-MMCT of ~220 to 370 μmol/kg (B/Ca-[B(OH) 4 −/DIC] sw ) and 280 to 430 μmol/kg (B/Ca-[B(OH) 4 −/HCO 3 −] sw) (Supplementary Table 5) is similar to change derived from a modeling study25, although absolute values are offset. Absolute Miocene DIC estimates here fall within the general estimates based on a range of approaches from modeling to geochemical reconstructions (Supplementary Fig. 12)35,90,91.

Regardless of the approach and location, our findings suggest the MCO DIC levels were higher than the early Miocene and late middle Miocene and DIC levels decreased across the MMCT by as much as 200–400 μmol/kg. Given the inherent uncertainties in these calculations, particularly in the [B] sw and the sensitivity of the B/Ca-[B(OH) 4 −/DIC] sw and [B(OH) 4 −/HCO 3 −] sw , the absolute values should be taken with caution. The trends recognized, however, because they are largely determined by B/Ca can be considered robust.

Although the Miocene B/Ca record from ODP site 761 presented in this study encompasses the short-term variations in δ13C (i.e., CM events), due to the lack of complimentary δ11B-pH records at the same resolution, we are unable to estimate relative DIC on these timescales, based on records from this study (Supplementary Fig. 3). Here we look to the T. trilobus B/Ca record from Malta which includes both B/Ca and δ11B-pH estimates across CM6 and shows a similar increase in B/Ca (Supplementary Fig. 3). Using the same approach as above and equations 4-7, we estimate relative change in DIC across CM6. A 17 μmol/mol increase in planktic B/Ca and corresponding increase of 0.08 in pH (Supplementary Fig. 3) across the CM6 event equates to a lowering of DIC by ~300 μmol/kg, although in absolute terms this remains poorly constrained at present. We do note that to fully quantify and assess surface DIC changes on these timescales, a suite of orbitally resolved records of the carbonate system (i.e., B/Ca and δ11B) across mid-Miocene carbon isotope excursions are needed.