A 15-million-year surface- and subsurface-integrated TEX 86 temperature record from the eastern equatorial Atlantic

. TEX 86 is a paleothermometer based on Thaumarcheotal glycerol dialkyl glycerol tetraether (GDGT) lipids and is one of the most frequently used proxies for sea-surface temperature (SST) in warmer-than-present climates. However, GDGTs are not exclusively produced in and ex-ported from the mixed layer, so sedimentary GDGTs may contain a depth-integrated signal that is also sensitive to local subsurface temperature variability. In addition, the correlation between TEX 86 and SST is not signiﬁcantly stronger than that to depth-integrated mixed-layer to subsurface temperatures. The calibration of TEX 86 to SST is therefore con-troversial. Here we assess the inﬂuence of subsurface temperature variability on TEX 86 using a downcore approach. We present a 15 Myr TEX 86 record from Ocean Drilling Program Site 959 in the Gulf of Guinea and use to elucidate the of the recorded TEX variability. Relatively high GDGT [ 2 / 3 ] ratio values from 13.6 that sedimentary GDGTs were deeper > 200 m) glacial–interglacial cyclicity, by benthic δ 18 O, while the variability within dinoﬂagellate assemblages and surface/thermocline temperature records (U k (cid:48) 37 and Mg / Ca) not by glacial–interglacial cyclicity. these observations are best explained by TEX 86 sensitivity to sub-thermocline temperature variability. We conclude that TEX 86 depth-integrated SST and


Introduction
Accurate tropical sea-surface temperature (SST) reconstructions are needed to assess oceanic heat distribution and ocean-atmosphere circulation in warmer-than-present climates, such as during the Pliocene. Tropical warming, even of a small magnitude, can lead to major changes in atmospheric circulation, with effects including intensification of monsoonal precipitation (Haywood et al., 2020; and warming of the extratropics via teleconnections (Barreiro et al., 2006). At present, the most frequently used proxies for past SSTs are Mg/Ca (foraminifer calcite), U k 37 (alkenones) and TEX 86 (glycerol dialkyl glycerol tetraethers; GDGTs). Each proxy has specific confounding factors and calibration issues. For example, Mg/Ca pa-leothermometry requires corrections for the Mg/Ca ratio of seawater and diagenesis (Dekens et al., 2002;Evans et al., 2016), and the U k 37 proxy is insensitive to temperatures above ∼ 28 • C (Müller et al., 1998). TEX 86 is potentially one of the best proxies to reconstruct SSTs above 28 • C, but one critical concern is that GDGTs are produced by marine Thaumarchaeota throughout the water column and often dominantly below the mixed layer (Ingalls et al., 2006;Kim et al., 2015;Lengger et al., 2014;Shah et al., 2008).
Although GDGTs may be dominantly produced at depths below the mixed layer, the efficiency of organic carbon export towards the sediment decreases exponentially with depth below the photic zone as shown with sediment trap data and simulations (e.g., Martin et al., 1987;Middelburg, 2019), which for GDGTs is supported by in situ measurements (e.g., Wüchter et al., 2005). However, as the production of GDGTs in the mixed layer is generally low relative to production within and below the thermocline, GDGTs produced below the mixed layer might still comprise a significant component of the sedimentary assemblage (e.g., Taylor et al., 2013;Ho and Laepple, 2017). The export depth of downcore GDGTs can potentially be traced with the fractional abundance of GDGT-2 to GDGT-3 (GDGT[2/3] ratio; Hernández-Sánchez et al., 2014;Pearson et al., 2016;Taylor et al., 2013;Villanueva et al., 2015). GDGT[2/3] values are typically low (< 5) in suspended particulate matter (SPM) sampled from shallow (< 100 m) water and sharply increase to > 25 in deeper water (Hurley et al., 2018). This division is currently best explained by depth differentiation of Thaumarcheotal ecotypes with different cyclization patterns (Kim et al., 2016;Villanueva et al., 2015). Core tops do not have GDGT[2/3] values indicative of a purely "deep" Thaumarcheotal GDGT source (Fig. 1). Instead, they show a gradual transition from low to moderately high values, which can be attributed to the integration of GDGTs from a range of shallow to intermediate depths. Likewise, TEX 86 (Eq. 1) may reflect a depthintegrated surface to subsurface signal (SubST) that is sensitive to temperature variability at and below the surface ocean (Kim et al., 2008(Kim et al., , 2012Schouten et al., 2002;Tierney and Tingley, 2015).
[GDGT-2 + GDGT-3 + cren ] [GDGT-1 + GDGT-2 + GDGT-3 + cren ] (1) Indeed, several multi-proxy paleostudies suggest that dissimilarities between TEX 86 -based SST reconstructions and other temperature records may be explained by a SubST signal recorded in TEX 86 variability (e.g., Lopes dos Santos et al., 2010;McClymont et al., 2012;Rommerskirchen et al., 2011;White and Ravelo, 2020). Together with the GDGT[2/3] distribution in core tops (Fig. 1), these downcore records support the use of a SubST calibration. However, the targeted depth interval and weight distribution of the temperature integration remain the subject of discussion Laepple, 2016, 2017;Tierney et al., 2014). The correlation of core top TEX 86 to SST is not significantly stronger than to SubST, as obtained from a wide range of depth intervals (Ho and Laepple, 2016;Schouten et al., 2002). Therefore, the sensitivity of sedimentary TEX 86 to temperature cannot be inferred from variability in the spatial dimension, but variability in the temporal dimension may provide solutions. Here, we present a new 0-15 Ma TEX 86 record from Ocean Drilling Program (ODP) Site 959 in the eastern equatorial Atlantic. We evaluate the TEX 86 and GDGT[2/3] evolution on million-year and glacial-interglacial timescales during the late Pliocene. We complement these records with U k 37 SST estimates and dinoflagellate cyst (dinocyst) assemblages and compare them to the Site 959 benthic δ 18 O record (Norris, 1998a;van der Weijst et al., 2020) and Mg/Ca-based SST and thermocline temperature records (van der Weijst et al., 2022) to elucidate the local source of downcore TEX 86 index variability and, consequently, paleoclimatic change.  Mascle et al., 1996). It is presently located below the eastward-flowing Guinea Current, which originates from the North Equatorial Counter Current and the Canary Current (Norris, 1998a). SST currently varies between 25.3 and 28.6 • C on a seasonal timescale , and surface salinity ranges between 34.6 and 35.0 . The water column is characterized by a shallow thermocline, with the 20 • C isotherm depth annually varying between ∼ 40 and 60 m in response to coastal upwelling van der Weijst et al., 2022). A minor upwelling event occurs in boreal winter, and a longer and stronger event occurs in boreal summer (Djakouré et al., 2017;Verstraete, 1992;Wiafe and Nyadjro, 2015).
We studied the 0-160 m interval of Hole C and 162-194 m interval of Hole A. The lithology gradually changes from nannofossil/foraminifer chalk at the bottom to nannofossil/foraminifer ooze at the top (Mascle et al., 1996). The age model is based on nannofossil biostratigraphy  and benthic δ 18 O stratigraphy (Wagner, 1998) between 0 and 23 meters below seafloor (m b.s.f.), benthic δ 18 O stratigraphy between 33 and 46 mbsf (van der Weijst et al., 2020), astronomical tuning of high-resolution X-ray fluorescence data (Vallé et al., 2016) between 51 and 100 mbsf, and a spline function through planktic foraminifer and nannofossil biostratigraphy between 162 and 194 mbsf (Norris, 1998b;Shafik et al., 1998). We plot the data between 21.35 and 97.43 mbsf on the revised meters composite depth (rMCD) scale of Vallé et al. (2016). Because outside this interval a splice is unavailable, we constructed a revised meters below sea floor (rMBSF) scale in coherence with Vallé The core was sampled at a relatively low resolution between 0 and 15 Ma for biomarker analysis and at a higher resolution in the late Pliocene (2.7-3.5 Ma) interval for multiproxy reconstructions. A total of 219 samples were extracted for TEX 86 , of which 105 are closely spaced in the 2.7-3.5 Ma interval. U k 37 data were generated for 79 of these samples, of which 52 are from the 2.7-3.5 Ma interval. Additionally, we performed palynological assessments on 48 late Pliocene splits of samples that were also used for biomarker analysis.

Biomarker proxies
At Utrecht University, sediment samples were freeze-dried, and the outer surfaces were removed using a clean knife to prevent contamination. The sediment was finely ground using a mortar and pestle. In the process, a split of coarsely ground material was taken for palynology. Biomarkers were extracted at 100 • C and 7.6×10 6 Pa using a dichloromethane (DCM) : methanol (9 : 1, v/v) solvent mixture in an accelerated solvent extractor (ASE 350, Dionex). The total lipid extract was separated over an activated Al 2 O 3 column into apolar, ketone and polar fractions using hexane : DCM (9 : 1, v/v), hexane : DCM (1 : 1, v/v) and DCM : methanol (1 : 1, v/v) as eluents, respectively. For GDGT analysis, a C46 GDGT (99 ng) standard was added to the polar fraction before it was filtered over a 0.45 µm PTFE filter and dissolved in hexane : isopropanol 99 : 1 (v/v). Measurements were performed on an ultrahigh-performance liquid chromatographmass spectrometer (UHPLC-MS) using the method of Hopmans et al. (2016). Samples with a branched to isoprenoid tetraether (BIT) index of > 0.3 could be influenced by terrestrially produced GDGTs (Hopmans et al., 2004;Weijers et al., 2006) and were excluded from further analysis (n = 3). All samples were within the reliable range of the methane index (MI; Zhang et al., 2011) and ring index (RI; , suggesting no appreciable inputs from methanogenic archaea or other non-thermal factors. We will present our results using various calibrations, including the conventional approach of calibrating TEX 86 to SST, but also using various models using the relation between TEX 86 and various depth integrations. The ketone fraction was dissolved in hexane and analyzed on a gas chromatograph (GC) coupled to a flame ionization detector (GC-FID, Hewlett Packard 6890 series) equipped with a CP-Sil 5 fused silica capillary column (25 m×0.32 mm; film thickness 0.12 µm) and a 0.53 mm precolumn. Samples were injected on-column at 70 • C with helium as a carrier gas and a flow rate of 2 mL min −1 . The oven program was as follows: 70 • C for 1 min, then ramped to 130 • C at 20 • C min −1 , then to 320 • C at 4 • C min −1 , and then held isothermal for 10 min. U k 37 values were calculated from the fractional abundances of the C 37:2 and C 37:3 alkenones following Prahl and Wakeham (1987) and calibrated to SST using the calibration of Müller et al. (1998). We chose not to adopt the Bayspline approach (Tierney and Tingley, 2018) or other tropical calibrations (e.g., Sonzogni et al., 1997, DSR II) as they may amplify noise into signal because of difficulty in accurately determining the very small C 37:3 peak areas at the high range of SST (Herbert et al., 2020).

Palynology
Coarsely crushed freeze-dried samples were spiked with a known amount of Lycopodium clavatum spores to quantify absolute palynomorph (Stockmarr, 1972) and treated with 30 % HCl (2×) and 40 % cold HF (2×) to remove carbonates and silicates. The residue was sieved over nylon mesh, and from the 10-250 µm fraction, the lighter organic fraction was separated from the heavy mineral fraction by suspension. The remaining palynological residue was mounted on glass microscopic slides with glycerine jelly. Samples were analyzed using an optical microscope under 400× magnification and approximately 300 dinoflagellate cysts (dinocysts) were counted per sample. Taxonomy follows that of Williams et al. (2017). Capisocysta was only identified at the genus level, because known Capisocysta species are best identified by the number of antapical plates (Head, 1998), which is complicated by the common dissociation of the hypocystal plates. We calculate the relative abundances of a taxon based on the total dinocyst sum. Dinocysts with affinities for cold-temperate waters in the modern ocean (Boessenkool et al., 2001;Zonneveld et al., 2013) were grouped as "cool species" and used for estimating surface water temperature. This species are Ataxodinium choane, Corrudinium devernaliae, Corrudinium harlandii, Impagidinium pallidum, Nematosphaeropsis labyrinthus and Pyxidinopsis reticulata. We also use the P /(P + G) ratio (Versteegh, 1994) as a paleoproductivity index. The P /(P +G) ratio quantifies the number of cysts produced by heterotrophic dinoflagellates, which are characterized by peridinioid-type tabulation (P ) over the total number of dinocysts, including those produced by autotrophic or mixotrophic species with a gonyaulacoid-type tabulation (G).

TEX 86 , U k 37 and GDGT[2/3] evolution at Site 959
We explore our TEX 86 results using several calibration approaches in Fig. 2, including an SST calibration using the exponential TEX H 86 (Kim et al., 2010). It should be noted that this calibration yields significant residuals in tropical regions and suffers from regression dilution (Tierney and Tingley, 2014). We also use a linear calibration model using a Bayesian method, a model specific for the tropical Atlantic and models assuming contributions from GDGTs produced below the mixed layer. Below, we first explore our TEX 86 record in the context of other GDGT-based indicators and other data using the TEX H 86 calibration for SST (Figs. 3-7), while subsurface calibrations are discussed in Figs. 7 and 8.
Between 15 and 11 Ma, TEX H 86 -based SST values fluctuate around an average value of ∼ 29 • C (Fig. 3), similar to the present-day local annual maximum surface temperature of 28.6 • C. A long-term gradual cooling starts around 11 Ma and is punctuated by two major cooling steps of ∼ 2 • C each, around 4.9 and 3.4 Ma, in which TEX H 86 drops well below the present-day local annual minimum surface temperature of 25.3 • C. The cooling trend seems to end in the late Pleistocene, where TEX H 86 records a ∼ 4 • C warming between 0.5 Ma and present, with a core top TEX H 86 of 24.0 • C. U k 37 is at saturation (∼ 28 • C) for most of the interval but registers some colder temperatures in the late Pliocene (27-28 • C) and drops below 27 • C after 1.8 Ma. The core top sample registers a temperature of 26.6 • C, which is within calibration error (1.5 • C; Müller et al., 1998) of the 27.5 • C modern mean annual SST . The GDGT[2/3] ratio is 7.0 in the core top. In the oldest part of the studied interval, GDGT[2/3] values vary between 4 and 5.5. Following an abrupt increase at ∼ 13.6 Ma, GDGT[2/3] varies largely between 5 and 8. After 1.7 Ma, GDGT[2/3] persistently increases to > 7.

Pliocene dinocyst assemblages
Palynological assemblages are dominated by dinocysts with only minor amounts of terrestrially derived pollen and spores. Dinocyst preservation is generally good. The assemblages are dominated by Brigantedinium spp. (19 %-60 %), closely followed by Spiniferites spp. (8 %-36 %) ( Fig. 4; see the Supplement for a complete overview of the assemblages). The interval between 3.33 and 3.16 Ma was sampled at a higher resolution (∼ 5 kyr on average) and displays strong variability of the negatively correlated Brigantedinium spp. and Spiniferites spp. on timescales < 10 kyr.
Other relatively abundant and consistently present taxa are Operculodinium spp., Pentapharsodinium dalei and Impagidinium spp. Peak abundances of Lingulodinium spp. were recorded at ∼ 3.26 Ma (Lingulodinium machaerophorum) and at 2.77 Ma (Lingulodinium hemicystum). Capisocysta sp. is consistently present in the lower part of the interval, with peak abundances (> 20 %) at 3.32 and 3.20 Ma, but has its last occurrence at 3.17 Ma. There is general coherence between TEX H 86 -SST variability and the relative abundance of the "cool species" group, which is dominated by N. labyrinthus. On average, the dinocysts with affinities for colder-temperate waters make up 6 % of the total assemblage, but peak abundances of 10 %-13 % occur during the M2 glacial. The P /(P + G) paleoproductivity index is 0.48 on average, with minimum and maximum values of 0.23 and 0.70, respectively (Fig. 4). P /(P + G) values are primarily driven by the abundance of Brigantedinium spp. The P /(P + G) ratio shows no significant relationship with SST as derived from TEX H 86 (R 2 = 0.05, p = 0.14 in linear regression analysis).

Identifying the source depth of the TEX 86 signal
The 15-million-year TEX H 86 record at Site 959 shows a cooling from > 30 • C between 11 and 15 Ma to < 20 • C between 0 and 3.3 Ma (Fig. 3). The U k 37 record, on the other hand, is near or at proxy saturation (> 27 • C) for most of the interval, with some lower values after 1.8 Ma. Assuming that both proxies reflect SST, the relatively small offset between U k 37 and TEX H 86 -SST offset in the core top sample (2.7 • C) could be explained by a combination of calibration errors  and seasonality. However, the average downcore offset between TEX H 86 and U k 37 is 5.4 • C and peaks at 7-8 • C during the late Pliocene interglacials (Figs. 5 and 6), which demonstrates that TEX 86 underestimates past SST at this site, assuming that U k 37 is a more accurate representation of SST. This could potentially be resolved by using an alternative regional SST calibration. The BAYSPAR-SST calibration (Tierney and Tingley, 2014) accounts for spatial variability of the TEX 86 -SST relationship and is the only available calibration that raises absolute SST estimates. However, similar to TEX H 86 (Kim et al., 2010) and low-latitude Atlantic (Zhang et al., 2018) calibrations, it likely overestimates SST variability (Fig. 2), indicating that the calibration slope is too steep. This is potentially a problem inherent to the calibration of TEX 86 to SST and could reflect downcore TEX 86 sensitivity to SubST instead of SST variability (Ho and Laepple, 2016). In the following sections, we explore multi-proxy data from Site 959 in search of indications for the water depth the downcore TEX 86 signal is reflecting.

GDGT[2/3] values
A substantial contribution of GDGTs from below the oceanic surface layer is supported by the relatively high GDGT[2/3] ratio, which fluctuates between 5 and 8 in the majority of the record (Fig. 3). Such values are rarely observed in modernday shallow water (< 100 m) SPM and surface sediments (e.g., Besseling et al., 2019;Hernández-Sánchez et al., 2014;Hurley et al., 2018). These values are best explained as a mixed signal of GDGTs from the upper 200 m and intermediate waters (Fig. 1c), i.e., likely dominantly from below the mixed layer (∼ 50 m) at Site 959. A cross-plot of TEX H 86 -SST and GDGT[2/3] (Fig. 7) shows that most data are clustered around a positive regression slope, but the oldest (> 13.6 Ma) and youngest (< 1.7 Ma) data plot outside this main cluster, potentially signaling systematic changes in GDGT production and export depth. Both 13.6 and 1.7 Ma mark the onset of a strong GDGT[2/3] increase (Fig. 3), signaling systematic deepening of the source of the sedimentary GDGT signal.
The ratio of GDGT[2] to GDGT[3] in samples does not directly affect the TEX 86 index because they are both included in the denominator and divisor (Eq. 1). If the abundance of GDGT[2] increases relative to that of GDGT[3] with export depth, then samples with a higher GDGT[2/3] ratio should be associated with lower TEX 86 values. The positive relationship between TEX 86 and GDGT[2/3] (Fig. 7) is therefore an interesting feature. Taylor et al. (2013) found a similar relation in the calibration dataset and in the paleodomain. Based on their analyses, they conclude that water depth is the dominant control of GDGT[2/3] ratios in sediments, recently supported by SPM analyses (Hurley et al., 2018). They also conclude that the positive correlation between SST and GDGT[2/3], notably the opposite relation to that expected based on homeoviscous adaptation, is an oceanographical artifact. Rather, the positive relation suggests that the long-term TEX 86 trends are not driven by GDGT export depth but reflect overall temperature changes instead. Alternatively, non-thermal factors such as water column oxygenation and nutrient supply may have influenced GDGT cyclization: high nutrient availability and ammonia oxidation rates have been linked to low TEX 86 values Park et al., 2018). However, following the GDGT[2/3] shift at 13.6 Ma, TEX H 86 does not systematically change until ∼ 11 Ma (Fig. 3). Moreover, late Pliocene dinocyst assemblages (Fig. 4) are characterized by a highly variable P /(P + G) ratio, indicating that upwelling and nutrient supply were highly variable on sub-Milankovitch timescales. In contrast to TEX H 86 (Fig. 5), the P /(P + G) ratio shows no glacial-interglacial variability. It is therefore unlikely that TEX 86 variability was primarily driven by non-thermal factors.

Surface, thermocline and sub-thermocline temperature variability
The late Pliocene TEX H 86 record follows a different evolution than SST records of U k 37 and the recent Mg/Ca of the surface-dwelling Globigerinoides ruber (Mg/Ca (G. ruber) ) and thermocline-dwelling Neogloboquadrina dutertrei (Mg/Ca (N. dutertrei) ) generated on the same samples (Fig. 5). The Mg/Ca records (van der Weijst et al., 2022) are based on very well-preserved foraminifera and apply a constant Pliocene correction for the Mg/Ca of seawater, a species-specific core top calibration with a constant partitioning coefficient and embedded dissolution correction (Dekens et al., 2002). Whereas TEX H 86 decreases by ∼ 2 • C, U k 37 and Mg/Ca-based SSTs show no significant net change in this interval. Moreover, although Mg/Ca and U k 37 register some cooling during the late Pliocene glacials, their variability seems not primarily driven by glacial-interglacial cyclicity, as is best observed in the high-resolution interval from M2 to KM2 (Fig. 5). This is also true for the relative abundance of dinocysts with affinities for colder temperatures (Figs. 4 and 5). Furthermore, around the KM2 and M2 glacials, TEX H 86 cooling leads to increasing abundances of cool dinocysts, as well as Mg/Ca (G. ruber) and U k 37 cooling by ∼ 5-10 ka (Fig. 6), which underscores the independent evolution of these records.
At sites where TEX 86 is suspected to be affected by GDGTs produced below the surface ocean, TEX H 86 has occasionally been interpreted to reflect thermocline temperature variability (e.g., Lopes dos Santos et al., 2010). Thermocline temperatures are lower and can be highly variable through time, which potentially explains both absolute SST underestimation and amplitude overestimation in specific TEX H 86 records. However, both could also be explained as a general artifact of an overestimated calibration slope that results from calibrating a depth-integrated GDGT signal to SST Laepple, 2015, 2016 (Figs. 3 and 5). Benthic δ 18 O is a faithful recorder of glacial-interglacial cyclicity because it registers a combined signal of deep ocean temperature and polar ice sheet volume. Could TEX 86 at Site 959 also be sensitive to glacial-interglacial cyclicity in deeper waters? At the onset of the M2 and KM2 glacials, the glacial expression of TEX H 86 leads δ 18 O by ∼ 5 kyr (Fig. 6), which indicates that, despite the similarities between TEX H 86 -SST and benthic δ 18 O, these proxies record variability in different water masses. Bottom waters at Site 959 are predominantly ventilated by North Atlantic Deep Water (NADW), whereas depths below the thermocline are occupied by South Atlantic Central Water (SACW; Fig. 8). SACW originates from sur-face waters in the South Atlantic gyre, which at depth mixes with Antarctic Intermediate Water (AAIW). It is therefore sensitive to climate variability in the mid-latitudes in the Southern Hemisphere and, like benthic δ 18 O, sensitive to higher-latitude climate change. Although our record samples SACW, AAIW dynamics have previously been suggested to be the main driver of TEX 86 variability in the Arabian Sea (Huguet et al., 2006) and in the southeast Atlantic (Rommerskirchen et al., 2011). If TEX 86 at Site 959 records a mixed signal from surface and subsurface waters, as is supported by the relatively high GDGT[2/3] ratio ( Figs. 1 and 3), it is likely sensitive to SACW dynamics, which explains the strong coherence with the benthic δ 18 O record. Collectively, these observations suggest that the downcore TEX 86 record at Site 959 is substantially affected by temperature variability below the surface ocean and that calibration to SubST is more appropriate than to SST.

Calibration of the depth-integrated TEX 86 record
Several TEX 86 calibrations for subsurface temperature (TEX 86 -SubST) are available, with each assuming, among other things, a different depth integration of the water column (Fig. 8). The BAYSPAR-SubST (Tierney and Tingley,  The oldest (> 13.6 Ma; dark blue) and youngest (< 1.7 Ma; dark red) samples largely stand out from the majority of the data. Linear regression lines determined from full dataset (stippled line) and samples dated between 13.6 and 1.7 Ma (solid line; 88 % samples, R 2 = 0.37, p < 0.001).

2014) and TEX H
86 -SubST (Kim et al., 2012) calibrations both target the upper 200 m of the water column, but the weight of the BAYSPAR-SubST calibration is centered at relatively shallow depths compared to the linearly weighted TEX H 86 -SubST calibration (Fig. 8). The HL16-SubST (Ho and Laepple, 2016) calibration targets the upper 950 m. The validity of this depth interval was questioned (Tierney et al., 2017), but with peak calibration weight at 100-350 m (Fig. 8), it is compatible with the documented vertical distribution of Thaumarcheotal cell counts and GDGT concentrations (e.g., Hernández-Sánchez et al., 2014;Schouten et al., 2012;Wuchter et al., 2005;Sintes et al., 2016), leaving export efficiency as the major uncertainty. Moreover, based on the large proportion of core tops with higher (> 5) GDGT[2/3] values in the global calibration set (Fig. 1), a slightly deeper calibration target (including depths below 200 m) is arguably preferable. This deeper calibration target lowers the slope of the calibration, which dampens reconstructed Miocene-Pliocene temperature variability at Site 959 (Fig. 2).
The ratio between temperature change in the surface and subsurface ocean is 1 : 1 when averaged across many sites and on longer timescales (Ho and Laepple, 2016). Under the assumption that the (integrated) depth of GDGT sourcing and the structure of the water column was stable, a depthintegrated TEX 86 record might still be of value to assess SST variability. Depending on the calibration, the magnitude of Late Miocene to Pleistocene cooling at Site 959 is 10 • C (BAYSPAR-SubST), 7 • C (TEX H 86 -SubST) or 5 • C (HL16- SubST; Fig. 8). Assuming a 1 : 1 ratio, the long-term SST trend should be of a similar magnitude. Currently available tropical SST estimates on the studied timescale are either compromised by saturation of the U k 37 proxy or based on TEX 86 and are therefore not suitable for independent comparison. However, global mean surface temperature estimates based on benthic δ 18 O data suggest a ∼ 5 • C cooling across the same interval (Hansen et al., 2013;Tierney et al., 2020). In the absence of major oceanographic changes, it is unlikely that the local cooling trend at Site 959 was larger, because temperature variability is generally lower in the tropics compared to high latitudes where deep-ocean waters derive from (i.e., polar amplification). This suggests that the BAYSPAR-SubST and TEX H 86 -SubST calibrations may overestimate the magnitude of long-term cooling (Fig. 8). From the currently available calibrations, HL16-SubST best approximates both expected glacial-interglacial variability and multi-million year trends (Figs. 2 and 8) and could potentially be used in combination with other SST records and/or modern water column data to correct for the local SST-SubST offset (Fig. 2). It is possible that the 1 : 1 ratio between SST and SubST did not persist across the entire record, e.g., due to changes in SACW production. To further improve tropical SST estimates in warmer climates, it should be further explored under which conditions a TEX 86 record that is driven by subsurface variability can be used to reconstruct SST variability and how to correct for the surface-subsurface temperature offset.

Late Neogene TEX 86 as an intermediate ocean signal at Site 959: new insights into M2 glacial inception
It is currently unclear what caused M2 glacial inception, but hypothesized mechanisms include declining atmospheric CO 2 levels (Berends et al., 2019;Dolan et al., 2015) and reduced latitudinal heat transport by the North Atlantic Current (De Schepper et al., 2009 in the Northern Hemisphere (NH) and by Indonesian throughflow (De Vleeschouwer et al., 2018) in the Southern Hemisphere (SH). Multi-proxy data from Site 999 in the Caribbean Sea show that declining CO 2 levels lagged the benthic δ 18 O glacial expression at the onset of M2 (de la Vega et al., 2020). This indicates that declining CO 2 levels did not initiate the M2 glacial, although they might have affected its intensity and duration. Because the TEX 86 record at Site 959 is affected by SACW, it can potentially be used to study the connection between our low-latitude site and higher latitudes. The TEX 86 and benthic δ 18 O records at Site 959 correspond closely but are out of sync during the M2 glacial (Fig. 6). Both records, generated on the same samples, display a double peak around M2, but TEX 86 leads δ 18 O by ∼ 5 kyr. Whereas the depthintegrated TEX 86 signal at Site 959 is sensitive to SH con- Figure 9. Exploring the chronology around the M2 glacial with records aligned on peak glacial benthic δ 18 O values (original age models in Fig. S2). (a) Daily insolation at 60 • N/S during NH/SH summer solstice (Laskar, 1990). (b) East equatorial Atlantic SST at Site 959 (this study) and Site 662 (Herbert et al., 2010) (Lisiecki and Raymo, 2005) and 999 (Haug et al., 2001;de la Vega et al., 2020), aligned to LR04 (Lisiecki and Raymo, 2005) (Berends et al., 2019).
We further explore the relative chronology around M2 in Fig. 9. The age models were aligned with peak glacial δ 18 O values, and the same records are plotted in the original age models in Fig. S1 in the Supplement. Although it is difficult to discriminate between signal and noise at this resolution, some patterns seem to emerge from this compilation.  Fig. 9, this intermediate ocean cooling occurred during a SH summer insolation minimum, suggesting that astronomical forcing may have amplified Southern Ocean cooling and tipped the system to glacial conditions. The sensitivity of the Antarctic ice sheet to astronomical forcing is also reflected by the 3.3 Ma onset of ∼ 20 kyr cyclicity in ice-rafted debris (IRD) at East Antarctic IODP Site U1361, with peak IRD mass accumulation rates at SH summer insolation minima (Patterson et al., 2014). Based on the chronology of Fig. 9, SH cooling led atmospheric CO 2 drawdown by ∼ 15-20 kyr. Lags of a similar magnitude were not unusual during Pleistocene glacial inceptions (Petit et al., 1999). Marginal Antarctic diatom assemblages and lithology suggest that considerable sea ice expansion occurred during M2, which could have inhibited CO 2 outgassing from the deep ocean (Ishino and Suto, 2020;McKay et al., 2012). Moreover, a slight increase in dust-mediated iron fertilization of the Southern Ocean may have facilitated export productivity (Martínez-Garcia et al., 2011). In combination with reduced mixing between Northern Component Water and Southern Component Water (van der Weijst et al., 2020), these processes would have led to increased carbon storage in the deep ocean in response to M2 cooling. Higher-resolution records and correlations would be required to fully test this hypothesis.

Conclusions
Several lines of evidence show that the TEX 86 record at Site 959 in the eastern equatorial Atlantic is best explained as a depth-integrated signal that is substantially affected by temperature variability below the thermocline. Relatively high GDGT[2/3] ratios since ∼ 13.6 Ma indicate that the sedimentary GDGTs were partly sourced from deeper (> 200 m) waters. Moreover, late Pliocene multi-proxy data show that TEX 86 variability is predominantly driven by glacial-interglacial variability, in contrast to other SST and thermocline records (U k 37 and Mg/Ca (G. ruber) SST and Mg/Ca (N. duterrei) ) and dinocyst abundances. The TEX 86 record strongly resembles the benthic δ 18 O record, indicating that it dominantly records temperature changes in a subthermocline water mass. At Site 959, intermediate depths are occupied by SACW, which derives from the mid-latitudes of the South Atlantic. According to the present-day water column composition, a substantial contribution of GDGTs from deeper (> 200) waters is needed to obtain TEX 86 sensitivity to SACW. This favors a subsurface calibration that integrates across a wider range of water depths, such as HL16-SubST (Ho and Laepple, 2016), in which calibration weight peaks 100-350 m. This calibration target interval is compatible with pelagic Thaumarcheotal cell counts and GDGT concentrations and with core top GDGT[2/3] values in the global calibration set.
Even if TEX 86 mainly reflects SubST variability, it may be used to reconstruct past SSTs if the temporal relationship between SST-SubST in a certain region is well-understood. Assuming a 1 : 1 relationship between long-term SubST and SST trends, our TEX 86 record suggests 5 • C tropical SST cooling between the Late Miocene and Pleistocene when calibrated with the Ho and Laepple (2016) subsurface calibration. Additionally, TEX 86 is also highly informative as a SubST proxy, because the sensitivity of TEX 86 to SACW at Site 959 offers a chance to explore connections between high-and low-latitude climate variability. The TEX 86 and benthic δ 18 O records at Site 959 are highly similar, but TEX 86 is sensitive to high-latitude SH conditions through AAIW, whereas the benthic δ 18 O record is primarily forced by high-latitude NH conditions through NADW. A ∼ 5 kyr lead of TEX 86 relative to δ 18 O at the onset of the M2 glacial stage indicates a lead of SH over NH cooling. Albeit limited by the resolution of the datasets, multi-proxy data from Site 959, Site 999 (Caribbean Sea) and Site 662 (eastern equatorial Atlantic) aligned based on LR04 peak glacial δ 18 O values indicate that SH cooling also led tropical ocean cooling and CO 2 levels by 15-20 kyr, suggesting that CO 2 drawdown was a consequence of glacial conditions. Instead, glacial expansion and SACW cooling at SH insolation minima suggest that orbital forcing played a pivotal role in M2 glacial inception. Data availability. New Site 959 data are available as the Supplement to this paper and have been submitted to the PANGAEA online data repository.
Author contributions. CMHvdW, FP, FS and AS designed the study. CMHvdW, KJvdL and TJTV generated the data. All authors contributed to interpretations. CMHvdW wrote the paper with input from FP, FS and AS and incorporated feedback from all authors.
Competing interests. At least one of the (co-)authors is a member of the editorial board of Climate of the Past. The peer-review process was guided by an independent editor, and the authors also have no other competing interests to declare.
Disclaimer. Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.