# Climate Evolution Through the Onset and Intensification of Northern Hemisphere Glaciation

E L Mcclymont, S L Ho, H L Ford, I Bailey, M A Berke, C T Bolton, S de Schepper, G R Grant, J Groeneveld, G N Inglis, et al.

► To cite this version:

E L Mcclymont, S L Ho, H L Ford, I Bailey, M A Berke, et al. Climate Evolution Through the Onset and Intensification of Northern Hemisphere Glaciation. Reviews of Geophysics, 2023, 61 (3), pp.e2022RG000793. 10.1029/2022rg000793. hal-04157006

HAL Id: hal-04157006

https://hal.science/hal-04157006

Submitted on 10 Jul 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Reviews of Geophysics®

## COMMISSIONED
MANUSCRIPT

10.1029/2022RG000793

†Deceased, 2022.

E. L. McClymont and S. L. Ho contributed equally to this work and should be considered joint first authors.

**Key Points:**

* The “stable” warm late Pliocene ~3.3–3.1 million years ago was a time of climate transition, especially in the southern hemisphere

* Ocean temperatures and ice sheets evolved asynchronously 3.3–2.4 Ma during the onset and intensification of Northern Hemisphere Glaciation

* Climate variability evolved in complex, non-uniform ways, most strongly expressed in northern mid-latitude sea-surface temperature records

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

E. L. McClymont and S. L. Ho,
erin.mcclymont@durham.ac.uk;
slingho@ntu.edu.tw

**Citation:**

McClymont, E. L., Ho, S. L., Ford, H. L., Bailey, I., Berke, M. A., Bolton, C. T., et al. (2023). Climate evolution through the onset and intensification of Northern Hemisphere Glaciation. *Reviews of Geophysics*, 61, e2022RG000793. https://doi.org/10.1029/2022RG000793

Received 5 OCT 2022
Accepted 26 APR 2023

**Author Contributions:**

**Conceptualization:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, S. De Schepper, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, B. Risebrobakken

© 2023. The Authors.

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

<figure><img src="image_2.png" /></figure>

## Climate Evolution Through the Onset and Intensification of Northern Hemisphere Glaciation

E. L. McClymont¹ Ⓓ ID, S. L. Ho² Ⓓ ID, H. L. Ford³ Ⓓ ID, I. Bailey⁴ Ⓓ ID, M. A. Berke⁵, C. T. Bolton⁶ Ⓓ ID, S. De Schepper⁷ Ⓓ ID, G. R. Grant⁸ Ⓓ ID, J. Groeneveld²,⁹ Ⓓ ID, G. N. Inglis¹⁰ Ⓓ ID, C. Karas¹¹ Ⓓ ID, M. O. Patterson¹² Ⓓ ID, G. E. A. Swann¹³ Ⓓ ID, K. Thirumalai¹⁴ Ⓓ ID, S. M. White¹⁵,† Ⓓ ID, M. Alonso-Garcia¹⁶ Ⓓ ID, P. Anand¹⁷ Ⓓ ID, B. A. A. Hoogakker¹⁸ Ⓓ ID, K. Littler¹⁹,²⁰ Ⓓ ID, B. F. Petrick²¹, B. Risebrobakken⁷ Ⓓ ID, J. T. Abell¹⁴ Ⓓ ID, A. J. Crocker¹⁰ Ⓓ ID, F. de Graaf³ Ⓓ ID, S. J. Feakins²² Ⓓ ID, J. C. Hargreaves²³, C. L. Jones¹⁰ Ⓓ ID, M. Markowska²⁴, A. S. Ratnayake²⁵ Ⓓ ID, C. Stepanek²⁶ Ⓓ ID, and D. Tangunan²⁷ Ⓓ ID

¹Department of Geography, Durham University, Durham, UK, ²Institute of Oceanography, National Taiwan University, Taipei, Taiwan, ³School of Geography, Queen Mary University of London, London, UK, ⁴Camborne School of Mines, University of Exeter, Penryn, UK, ⁵Department of Civil and Environmental Engineering and Earth Sciences, University of Notre Dame, Notre Dame, IN, USA, ⁶Aix Marseille University, CNRS, IRD, INRAE, CEREGE, Aix-en-Provence, France, ⁷NORCE Climate and Environment, NORCE Norwegian Research Centre, Bjerknes Centre for Climate Research, Bergen, Norway, ⁸Department of Surface Geosciences, GNS Science, Lower Hutt, New Zealand, ⁹Department of Geology, Hamburg University, Hamburg, Germany, ¹⁰School of Ocean and Earth Science, National Oceanography Centre Southampton, University of Southampton, Southampton, UK, ¹¹Universidad de Santiago de Chile, Santiago, Chile, ¹²Department of Geological Sciences and Environmental Studies, Binghamton University, Binghamton, NY, USA, ¹³School of Geography, University of Nottingham, Nottingham, UK, ¹⁴Department of Geosciences, The University of Arizona, Tucson, AZ, USA, ¹⁵Department of Earth and Planetary Sciences, University of California, Santa Cruz, CA, USA, ¹⁶Department of Geology, Universidad de Salamanca, Salamanca, Spain, ¹⁷School of Environment Earth and Ecosystem Sciences, The Open University, Milton Keynes, UK, ¹⁸The Lyell Centre, Heriot-Watt University, Edinburgh, UK, ¹⁹Department of Earth and Environmental Sciences, University of Exeter, Penryn, UK, ²⁰Environment and Sustainability Institute, University of Exeter, Penryn, UK, ²¹Institute of Geosciences, Kiel University, Kiel, Germany, ²²Department of Earth Science, University of Southern California, Los Angeles, CA, USA, ²³Blue Skies Research Ltd, Settle, UK, ²⁴Climate Geochemistry, Max Planck Institute for Chemistry, Mainz, Germany, ²⁵Department of Applied Earth Sciences, Uva Wellassa University, Badulla, Sri Lanka, ²⁶Department for Paleoclimate Dynamics, Alfred Wegener Institute—Helmholtz Centre for Polar and Marine Research, Bremerhaven, Germany, ²⁷School of Earth and Environmental Sciences, Cardiff University, Cardiff, UK

**Abstract** The Pliocene Epoch (~5.3–2.6 million years ago, Ma) was characterized by a warmer than present climate with smaller Northern Hemisphere ice sheets, and offers an example of a climate system in long-term equilibrium with current or predicted near-future atmospheric CO₂ concentrations (pCO₂). A long-term trend of ice-sheet expansion led to more pronounced glacial (cold) stages by the end of the Pliocene (~2.6 Ma), known as the “intensification of Northern Hemisphere Glaciation” (iNHG). We assessed the spatial and temporal variability of ocean temperatures and ice-volume indicators through the late Pliocene and early Pleistocene (from 3.3 to 2.4 Ma) to determine the character of this climate transition. We identified asynchronous shifts in long-term means and the pacing and amplitude of shorter-term climate variability, between regions and between climate proxies. Early changes in Antarctic glaciation and Southern Hemisphere ocean properties occurred even during the mid-Piacenzian warm period (~3.264–3.025 Ma) which has been used as an analog for future warming. Increased climate variability subsequently developed alongside signatures of larger Northern Hemisphere ice sheets (iNHG). Yet, some regions of the ocean felt no impact of iNHG, particularly in lower latitudes. Our analysis has demonstrated the complex, non-uniform and globally asynchronous nature of climate changes associated with the iNHG. Shifting ocean gateways and ocean circulation changes may have pre-conditioned the later evolution of ice sheets with falling atmospheric pCO₂. Further development of high-resolution, multi-proxy reconstructions of climate is required so that the full potential of the rich and detailed geological records can be realized.

**Plain Language Summary** Warm climates of the geological past provide windows into future environmental responses to elevated atmospheric CO₂ concentrations, and past climate transitions identify important or sensitive regions and processes. We assessed the patterns of average ocean temperatures and indicators of ice sheet size over hundreds of thousands of years, and compared to shorter-term variability (tens of thousands of years) during a recent transition from late Pliocene warmth (when CO₂ was similar to present) to the onset of the large and repeated advances of northern hemisphere ice sheets referred to as the “ice ages.”
**Data curation:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, C. T. Bolton, S. De Schepper, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, S. M. White, M. Alonso-Garcia, P. Anand, K. Littler, B. F. Petrick, B. Risebrobakken, J. T. Abell, A. J. Crocker, F. de Graaf, S. J. Feakins, J. C. Hargreaves, C. L. Jones, M. Markowska, A. S. Ratnayake, C. Stepanek, D. Tangunan

**Formal analysis:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, P. Anand, B. A. A. Hoogakker, K. Littler, B. F. Petrick, B. Risebrobakken

**Investigation:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, P. Anand, B. A. A. Hoogakker, K. Littler, B. F. Petrick, B. Risebrobakken, J. T. Abell, A. J. Crocker, S. J. Feakins, M. Markowska, C. Stepanek, D. Tangunan

**Methodology:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, C. T. Bolton, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, B. A. A. Hoogakker, K. Littler

**Project Administration:** E. L.

McClymont, S. L. Ho, H. L. Ford

**Resources:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, C. T. Bolton, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, P. Anand, B. A. A. Hoogakker, B. F. Petrick, B. Risebrobakken, J. T. Abell, A. J. Crocker, F. de Graaf, S. J. Feakins, J. C. Hargreaves, C. L. Jones, M. Markowska, A. S. Ratnayake, C. Stepanek, D. Tangunan

**Software:** S. L. Ho, G. R. Grant, K. Thirumalai

**Validation:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, C. T. Bolton, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, B. A. A. Hoogakker, K. Littler, B. Risebrobakken

**Visualization:** E. L. McClymont, S. L. Ho, G. R. Grant, G. E. A. Swann, K. Thirumalai, M. Alonso-Garcia, B. A. A. Hoogakker

**Writing – original draft:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, P. Anand, B. A. A. Hoogakker, K. Littler, B. Risebrobakken

We show that different regions of the climate system changed at different times, with some changing before the ice sheets expanded. The development of larger ice sheets in the Northern Hemisphere then impacted ocean temperatures and circulation, but there were many regions where no impacts were felt. Our analysis highlights regional differences in the timing and amplitudes of change within a globally-significant climate transition as well as in response to the current atmospheric CO₂ concentrations in our climate system.

# 1. Introduction

Over the last ~50 million years (Myr), there has been a long-term cooling trend in the Earth's climate, culminating in the transition to an “icehouse” climate state of bipolar and high amplitude glaciations which began in the Pliocene Epoch (5.3–2.58 million years ago, Ma) and was fully established across the Pliocene-Pleistocene transition (~2.58 Ma) (Figure 1) (Westerhold et al., 2020). The Late Pliocene thus offers us an important opportunity to study the characteristics and pacing of major global climate change, and allows us to consider whether regions or parts of the climate system are vulnerable or resilient to climate forcings and feedbacks.

Late Pliocene (3.6–2.586 Ma) climate was characterized by a sustained period (~200,000 years) of global warmth, the mid-Piacenzian warm period (mPWP, ~3.264–3.025 Ma). The mPWP has been an important window for collecting information about Earth's climate response to atmospheric carbon dioxide (CO₂) concentrations comparable to today and our near-future (~400 ppmv; Figure 1) (Gulev et al., 2021). The Late Pliocene is also marked by the development of more intense cold stages (glacials) and an increase in ice-rafted debris (IRD) reaching the North Atlantic and North Pacific Oceans (Figure 1) (Blake-Mizen et al., 2019; Jansen et al., 2000; Kleiven et al., 2002; Mudelsee & Raymo, 2005). These changes mark the final step in the overall trend toward a bipolar-icehouse climate state (Bailey et al., 2013; Westerhold et al., 2020), and have been variously referred to as the “onset” or the “intensification of Northern Hemisphere Glaciation” (oNHG or iNHG). However, the terms “oNHG” and “iNHG” have been inconsistently used in the literature. A wide window of time in which a gradual expansion of the ice sheets developed, starting as early as ~3.6 Ma (e.g. Mudelsee & Raymo, 2005), has been referred to as iNHG (e.g. Bolton et al., 2010, 2018; Kleiven et al., 2002). In contrast, iNHG has also described a narrower time window aligned with the stepped increase in IRD recorded in the high northern latitudes at ~2.7 Ma (e.g. Bailey et al., 2011; Maslin et al., 1998; Ravelo et al., 2004; Teschner et al., 2016). Others have referred to the changes at ~2.7 Ma as the “onset of major Northern Hemisphere Glaciation” (oNHG), to recognize that whilst ice was present and developing during the Pliocene (cf. Mudelsee & Raymo, 2005) the development of pronounced glacial cycles and larger continental ice sheets was a later significant climate shift (e.g. Bailey et al., 2013; Haug et al., 2005). Here, we make a distinction between (a) the longer-term transition toward bipolar glaciations, potentially beginning as early as ~3.6 Ma and including the mPWP (which we refer to as the onset of NHG, oNHG), and (b) the relatively abrupt and later expansion of Northern Hemisphere ice sheets which occurs from ~2.7 Ma (which we refer to as the intensification of NHG, iNHG).

One hypothesis to explain the iNHG is that the climate system crossed a threshold whereby a fall in CO₂ allowed large-scale glaciation in the Northern Hemisphere (Figure 1) (De Conto et al., 2008; Lunt et al., 2008). However, both the oNHG and iNHG occurred in the context of long-term shifts to meridional and zonal temperature gradients (e.g. Fedorov et al., 2015; Kaboth-Bahr & Mudelsee, 2022) as well as tectonic changes which can influence heat and moisture transport through the climate system, with potential impacts on ice-sheet growth (e.g. De Vleeschouwer et al., 2022; Karas et al., 2009; Knies, Mattingsdal, et al., 2014; Sánchez-Montes et al., 2020). It is important to determine the causes of the oNHG and iNHG, because they have been shown to affect the global climate system response to changes in the Earth's orbit and axial tilt (“orbital forcing”; Meyers & Hinnov, 2010; Turner, 2014; Westerhold et al., 2020). Orbital forcing impacts the distribution of energy received by the Earth from the sun through time and space, and has been shown to pace climate over tens and hundreds of thousands of years (e.g. Hays et al., 1976). A growing signal of global ice volume variability which could be linked to orbital forcing has been proposed to represent the development of a “deterministic” climate system associated with iNHG, whereby orbital forcing plays an important role in determining the climate system change (Meyers & Hinnov, 2010). In contrast, before iNHG the climate system is argued to be less strongly linked to orbital forcing, with a more stochastic (random) nature (Meyers & Hinnov, 2010). However, judging climate system response to forcing using only globally-integrated indicators, such as global ice volume, limits our assessment of the regional or local processes and feedbacks which may be facilitating the changes to external forcing as well as the underlying growth of the continental ice sheets.
**Writing – review & editing:** E. L. McClymont, S. L. Ho, H. L. Ford, I. Bailey, M. A. Berke, C. T. Bolton, S. De Schepper, G. R. Grant, J. Groeneveld, G. N. Inglis, C. Karas, M. O. Patterson, G. E. A. Swann, K. Thirumalai, S. M. White, M. Alonso-Garcia, P. Anand, B. A. A. Hoogakker, K. Littler, B. F. Petrick, B. Risebrobakken, J. T. Abell, A. J. Crocker, F. de Graaf, S. J. Feakins, J. C. Hargreaves, C. L. Jones, M. Markowska, A. S. Ratnayake, C. Stepanek, D. Tangunan

The Past Global Changes (PAGES) working group “PlioVAR” (Pliocene climate variability on glacial-interglacial timescales) (McClymont et al., 2017, 2020) has sought to compile and assess regional climate records spanning the Late Pliocene and the transition to the Early Pleistocene (i.e., the oNHG), to include both the mPWP and iNHG intervals (3.3–2.4 Ma). Our focus here is on records which have been recovered from marine sediment cores, including changes in sea-surface temperatures (SSTs), global ice volume and sea level. As SSTs can be reconstructed using a variety of approaches, our first research question was to consider whether proxy choice influenced our understanding of the climate signals we have examined:

(Q1) Are there proxy-specific differences in the records of SST change through time? Several proxy (indirect) methods of reconstructing SSTs are available, including assemblages and chemistry of both organic and inorganic remains of marine biota found in sediments (e.g., Dowsett et al., 2012; McClymont et al., 2020). Each proxy has different biological and environmental controls over how the temperature

<figure><img src="image_3.png" /></figure>

**Figure 1.** Key climate parameters across the Pliocene-Pleistocene transition. The uppermost panel delineates the transition from the Pliocene Epoch (3.6–2.58 Ma) to the Pleistocene Epoch (2.58–0.01 Ma) and the three time intervals used for our statistical analysis (P-1 3.3–3.0 Ma, P-2 3.0–2.7 Ma, and P-3 2.7–2.4 Ma; see Section 2.3). The time intervals defined as the onset of Northern Hemisphere Glaciation (oNHG), mid-Piacenzian warm period (mPWP) and the intensification of major Northern Hemisphere Glaciation (iNHG) are marked on panel (b), alongside the M2 glacial stage which has been suggested as a possible harbinger of more extensive glaciation (Westerhold et al., 2020). (a) Reconstructed atmospheric $pCO_2$ concentrations (de la Vega et al., 2020) including the ~280 ppmv threshold for extensive glaciation proposed by De Conto et al. (2008) which aligns with Pre-Industrial $pCO_2$ recorded in ice cores (Loulerge et al., 2008); (b) a global stack of benthic oxygen isotope ratios (Ahn et al., 2017) which reflect a combined signal of the deep-water temperature and ice volume (see Technical Box). Key marine isotope stages referred to in the text are marked, blue for the cold glacial stages, red for the warm interglacial stage KM5c; (c) North Atlantic ice-rafted debris (IRD) inputs at Sites U1307 (southern Greenland margin) and U1313 (southern North Atlantic) (Blake-Mizen et al., 2019; Lang et al., 2014), which may be used to indicate when ice sheets expanded and reached the sea. (d) Antarctic IRD inputs from Site U1361 (Patterson et al., 2014). MAR refers to the Mass Accumulation Rate of the sediment.
signal is recorded, including whether the signal is generated by photosynthesizers or grazers, the preferred water depth or season of production, and the preservation of the SST signal during transport through the water column and into the seafloor sediments. By understanding both the similarities and differences between SST proxy records, we may identify the processes influencing our temperature signals and gain a more detailed view of environmental change associated with the oNHG and iNHG.

The information from Q1 is then combined with evidence for changes in global ice volume, including IRD, and other records of climate change to address three key questions linked to climate forcing and responses associated with oNHG and iNHG:

(Q2) What were the characteristics of the mPWP interval? The mPWP is the most recent interval of geological time where pCO₂ and global temperatures were sustained above Pre-Industrial (by ~50–110 ppmv for pCO₂, by ~2.5–4.0°C for temperature) (de la Vega et al., 2020; Gulev et al., 2021). For this reason, the mPWP has long been a target for data synthesis and data-model comparison efforts (e.g., Dowsett et al., 2012, 2016; Haywood et al., 2020; McClymont et al., 2020; Salzmann et al., 2013), because it offers an opportunity to compare a warmer-than-modern climate with both modern observations and near-future projections (Burke et al., 2018; Gulev et al., 2021; Tierney et al., 2020). However, climate models of two interglacials within the mPWP show the potential for pronounced regional and seasonal sensitivity to variations to Earth's orbital configurations (Prescott et al., 2014), highlighting the importance of considering regional expressions of mPWP climates to better understand forcings and feedbacks on orbital timescales.

(Q3) What were the characteristics of the iNHG? The mPWP and iNHG mark a shift in the patterns of climate change over the last ~50 Ma, including a stronger global influence of complex high-latitude climate dynamics (Meyers & Hinnov, 2010; Turner, 2014; Westerhold et al., 2020). However, these analyses are focused on the integrated records of global ice volume and deep water temperatures recovered from stacks of benthic foraminifera stable oxygen isotope ratios (δ¹⁸O_benthic; Figure 1). Assessing a combination of globally-distributed surface ocean temperature and individual δ¹⁸O_benthic records offers opportunities to identify and explain long-term trends in mean climate state as well as climate variability (see Question 4) for specific regions and/or circulation systems.

(Q4) Did the amplitude or pacing of climate variability change with the iNHG? The climate system during the Pliocene and early Pleistocene was regulated on orbital timescales (10³–10⁵ kyr) (Lisiecki & Raymo, 2005; Shackleton et al., 1984). However, the dominance of orbital obliquity periods (linked to Earth's axial tilt, ~41 kyr) and the absence of a strong orbital precession signal (~19–23 kyr) is puzzling, as daily summer insolation intensity has strongly influenced changes in polar ice volume for the last ~800 kyr and is paced by precession (Hays et al., 1976). By exploring regional variability in ocean temperature and individual δ¹⁸O_benthic records we can explore whether global ice volume is key to the pacing of orbital-scale climate variability, and how that might have changed with the iNHG.

Here, we draw on the high temporal resolution, multi-proxy marine sediment data synthesis which was generated and evaluated by the PAGES-PlioVAR working group (Figure 2). In previous work we have focused on a single interglacial within the mPWP (KM5c, ~3.205 Ma; Figure 1) where well-constrained pCO₂ reconstructions allowed a detailed examination of the ocean temperature response to pCO₂ (McClymont et al., 2020). Here we apply the same stratigraphic and data evaluation protocols to the multi-proxy data over a longer time interval, to address the four questions outlined above.

We recognize that aspects of our data evaluation and analyses warrant technical explanations and justifications for a broad geophysics audience. In Section 2 we therefore first outline our approaches to ensuring sites meet our high-quality age control (Section 2.1) and briefly introduce both the climate records we have used (Section 2.2) and the statistical methods we employed to analyze our data compilation (Section 2.3). We direct readers who are interested in the details of the principles which underpin the climate reconstruction techniques we used to the Technical Box. The climate proxies we employ are indirect measures of climate variables, and through multi-proxy analysis we endeavor to tease out the most important signatures of past climate change. In particular, we employ both organic and inorganic proxies of ocean temperature change, which have a range of biological and environmental controls over their signatures: in Section 3 we address Q1 outlined above, by using three regions where we have the highest density of multi-proxy data to explore the differences and similarities expressed by these different approaches. In Section 4 we present our results examining changes in ocean temperature, ice volume and sea level through the late Pliocene to early Pleistocene (3.3–2.4 Ma). We first assess long-term
<figure><img src="image_4.png" /><figcaption>**Figure 2.** The PlioVAR-synthesized records of climate change spanning the 3.3–2.4 Ma interval, superimposed upon present-day mean annual sea-surface temperatures (SSTs). Four climate proxies for ocean temperature were assessed: $U'_{37}$, TEX$_{86}$, Mg/Ca in planktonic foraminifera, and $\delta^{18}$O in planktonic foraminifera. Benthic $\delta^{18}$O records were analyzed to assess changes in deep water temperature and global ice volume. The three regions which are discussed when comparing SST proxies (Section 3) are annotated, noting that upwelling sites are found in the Atlantic, Indian and Pacific Oceans. Not shown are the three sea-level records, since they assess global stacks of available $\delta^{18}$O$_{benthic}$. Sites which did not meet the PlioVAR thresholds for age control or temporal resolution of data are indicated in gray.</figcaption></figure>

changes (Section 4.1) and then examine glacial-interglacial variability (Section 4.2) using the “PlioVAR data sets” which met our data quality thresholds as outlined in Section 2.1. We then address Q2 (Section 4.3), Q3 (Section 4.4), and Q4 (Section 4.5), by synthesizing our findings in the context of published information on changes to ice sheets and regional climates.

Technical Box Reconstructing past ocean temperatures, ice volume and sea
level using marine sediment records

1. $U'_{37}$ temperature proxy: principles and interpretations.

Most of the SST records (17 sites) are based on the $U'_{37}$ proxy (Table 1), which is a ratio describing the relative distribution of specific organic molecules (the C$_{37}$ alkenones) synthesized by selected haptophyte algae (Marlowe et al., 1984; Prahl & Wakeham, 1987). Many Pliocene studies have calibrated $U'_{37}$ to mean annual SST using the linear calibrations of seafloor surface sediments (core-tops) (Müller et al., 1998) or laboratory cultures (Prahl & Wakeham, 1987), which are almost indistinguishable. However, for some regions, these calibrations do not align with recorded mean annual SST. Non-linear $U'_{37}$-temperature relationships are observed for high (i.e., warmest) $U'_{37}$ values (Conte et al., 2006; Pelejero & Calvo, 2003; Tierney & Tingley, 2018). Seasonality also likely impacts $U'_{37}$ reconstructions in the high latitudes of the North Pacific and North Atlantic Oceans (Tierney & Tingley, 2018). In our previous work on the interglacial KM5c (Figure 1), we showed that there was data-model agreement for both seasonal and mean annual SST reconstructions in the high latitudes, regardless of the calibration (McClymont et al., 2020). In the low-latitudes, the application of the non-linear BAYSPLINE calibration (Tierney & Tingley, 2018) elevated tropical SSTs relative to the Müller et al. (1998) core-top calibration (McClymont et al., 2020). Here, we present $U'_{37}$-temperatures calibrated using BAYSPLINE (Tierney & Tingley, 2018) to acknowledge the non-linearity of the $U'_{37}$-SST relationship at high values which is not captured by the Müller et al. (1998) calibration. However, we note that in doing so, there is increased uncertainty on the reconstructed SSTs above 24°C (from ~1.5 to ~4.4°C; Tierney & Tingley, 2018).
## 2. TEX<sub>86</sub> temperature proxy: principles and interpretations.

Only one published record (Site U1463) had TEX<sub>86</sub> temperature data meeting our age resolution criteria for analysis (Figure 1, Table 1). However, we draw on lower resolution TEX<sub>86</sub> records as secondary evidence in evaluating potential controls over other ocean temperature proxies for 4 other sites (Section 3).

The TEX<sub>86</sub> index is based on the distribution of isoprenoidal glycerol dialkyl glycerol tetraethers (isoG-DGTs) in marine sediments (Schouten et al., 2002), and has been correlated to surface or subsurface temperatures, where subsurface can be tens or a few hundred meters below the sea surface. Both linear (Schouten et al., 2002; Tierney & Tingley, 2014, 2015) or non-linear (Kim et al., 2010) calibrations have been used (see Inglis & Tierney, 2020, for an extensive review). Both calibration approaches yield similar values within the temperature range ~5–30°C. TEX<sub>86</sub> has been particularly useful where U<sup>K'</sup><sub>37</sub> values have reached their upper limit in warm waters (i.e., U<sup>K'</sup><sub>37</sub> = 1) (L. Li et al., 2011; O'Brien et al., 2014; van der Weijst et al., 2022; Y. G. Zhang et al., 2014). As for U<sup>K'</sup><sub>37</sub>, knowledge of past seawater chemistry is not required to calculate temperature values. However, care is required to ensure that data are not biased by contributions from archaea other than marine *Thaumarchaeota*, including terrestrial (Weijers et al., 2006), methanogenic (Inglis et al., 2015) and methanotrophic inputs (Y. G. Zhang et al., 2011). It has also been shown that changes to thermocline depth can influence TEX<sub>86</sub> values (Ford et al., 2015) (Section 3). We have screened all of the TEX<sub>86</sub> data included here using established indices for non-*Thaumarchaeota* inputs (e.g. BIT index, Ring index, Methane Index; Hopmans et al., 2004; Y. G. Zhang et al., 2011, 2016). We also evaluated the potential contribution from “deep-water” *Thaumarchaeota* (those living below the permanent pycnocline) using the GDGT-2/GDGT-3 ratios (following Rattanasriampaipong et al., 2022; K. W. Taylor et al., 2013). We then applied a spatially-varying linear Bayesian regression model (BAYSPAR, Tierney & Tingley, 2014) to calculate TEX<sub>86</sub>-temperature estimates assuming a surface ocean origin (0 m water depth) for the isoGDGTs.

## 3. Foraminifera Mg/Ca temperature proxy: principles and interpretations.

The Mg/Ca-temperature proxy is underpinned by the observation that Ca is preferentially substituted by Mg in the calcite crystal lattice of foraminifera shells at higher temperatures (e.g. Anand et al., 2003; Dekens et al., 2002; Lea et al., 1999; Nürnberg et al., 2000). The Mg/Ca-temperature relationship is exponential, such that the proxy is more sensitive at higher temperatures. Planktonic foraminifera species commonly used for surface water temperature reconstruction include *Globigerina bulloides*, *Globigerinoides ruber*, and *Trilobatus sacculifer* (renamed from *Globigerinoides sacculifer*, Spezzaferri et al., 2015). However, depending on the oceanographic setting, there may be species-specific preferences in the preferred water depth of their habitat or the seasonality of their productivity, which can influence the recorded temperature signal. For example, the tropical and subtropical species *G. ruber* and *T. sacculifer* have photosynthetic symbionts and are confined to the photic zone, with *T. sacculifer* occupying a slightly deeper depth habitat than *G. ruber* (e.g. Bé, 1980; Fairbanks et al., 1982; Rebotim et al., 2017). Species-specific Mg/Ca-temperature calibrations also highlight species-specific differences in the substitution of Mg (Anand et al., 2003; Elderfield & Ganssen, 2000) and its distribution within foraminifera tests (Anand & Elderfield, 2005; Brown & Elderfield, 1996; Spero et al., 2015). Calcification temperature is the dominant control on planktic foraminiferal Mg/Ca (e.g. Anand et al., 2003; Dekens et al., 2002; Lea et al., 1999; Nürnberg et al., 2000), despite secondary controls over Mg incorporation into foraminifera tests. These secondary controls may be corrected as part of Mg/Ca-temperature calculations, and include (a) salinity, (b) surface water pH, (c) the Mg/Ca ratio of seawater, and (d) partial dissolution of calcite at the seafloor by sediment porewaters (Dekens et al., 2002; Evans et al., 2016; Regenberg et al., 2006; Rosenthal et al., 2022).

Six foraminifera Mg/Ca data sets met our criteria for analysis, from *G. bulloides*, *G. ruber*, and *T. sacculifer* (Table 1). The PlioVAR working group has previously presented the originally-published SSTs from these data sets for the KM5c interval (McClymont et al., 2020). The KM5c synthesis thus included temperature values generated across a range of calibrations and corrections, which were then compared to new PlioVAR calculations using the Bayesian regression model “BAYMAG” (Tierney et al., 2019).
Here, we calculate ocean temperature from Mg/Ca data sets where continuous data had been generated from a single species across our target interval (3.3–2.4 Ma); we use an independent approach that maximized the fit between available core-tops and modern ocean temperatures (see Section 3 for results), while also maintaining consistency in our approaches to considering the secondary impacts on Mg/Ca incorporation noted in the previous paragraph, between sites and among species. We first accounted for the potential impact of dissolution at sites where modern-day carbonate saturation state is lower than 21.3 µmol kg⁻¹ (identified by Regenberg et al. (2014) as the critical value below which dissolution starts). To do this, we identified the modern-day saturation state using values stated in the original publications (Table 1), or using World Ocean Atlas data and the CO2calc software (Robbins et al., 2010). Second, where reductive cleaning was included in the preparation method, 10% was added to the original Mg/Ca values before calibration, to enable comparison with non-reductively cleaned samples (Barker et al., 2003; Khider et al., 2015; Rosenthal et al., 2004); Sites ODP 806 and 1143). Next, we converted all corrected Mg/Ca data to temperature using the Elderfield and Ganssen (2000) calibration (G. bulloides, Site U1313) and the Anand et al. (2003) calibrations (all other sites/species), which yielded the best fit to core-top values (Section 3).

Finally, we account for the effect of changes in Mg/Ca<sub>seawater</sub> on foraminiferal Mg/Ca. Shifts in Mg/Ca<sub>seawater</sub> operate over million-year timescales, and thus should not affect the relative variability which is our focus here. It has nevertheless been shown that there can be a non-linear impact of Mg/Ca<sub>seawater</sub> on the calculated absolute SSTs (and thus the ranges) which would impact our assessment of glacial-interglacial variability (Evans et al., 2016). Although the magnitude and method of Mg/Ca<sub>seawater</sub> correction remains uncertain (e.g. Rosenthal et al., 2022; White & Ravelo, 2020), comparison of Mg/Ca-based and clumped isotope-based temperatures in the Pliocene show that a modest correction is best supported (Meinicke et al., 2021). Furthermore, since the application of a variety of Mg/Ca-based calibrations gave temperature estimates consistent within the calibrations' uncertainty (<±1°C) (Rosenthal et al., 2022), the Mg/Ca<sub>seawater</sub> calibration choice has minimal impact on the timescales we are targeting here. We therefore used the Evans et al. (2016) Mg/Ca<sub>seawater</sub>-sensitive temperature calibration, in conjunction with the BAYMAG seawater reconstruction (Tierney et al., 2019) to determine the magnitude of the Mg/Ca<sub>seawater</sub> correction which we then applied to our records. *T. sacculifer* Mg/Ca values were corrected upward by 10.3% to align with *G. ruber* before calculating the Mg/Ca<sub>seawater</sub> correction, as per Evans et al. (2016). The Mg/Ca<sub>seawater</sub> correction (about +0.5°C) was then added to each calculated sample temperature. Overall, our approach yields more reasonable core-top SSTs than when only using the Evans et al. (2016) calibration, the latter giving unrealistically cold core-top temperatures if interpreted as reflecting SST. A summary of the original published Mg/Ca data and the impact of our PlioVAR corrections is shown in Figure S2 in Supporting Information S1.

### 4. Foraminifera δ¹⁸O records for surface and deep ocean properties.

Seven of the PlioVAR synthesis sites provided planktonic foraminiferal δ¹⁸O records (δ¹⁸O<sub>planktonic</sub>; Table 1). These δ¹⁸O<sub>planktonic</sub> records provide information about surface and near-surface temperatures and the δ¹⁸O composition of the seawater in which those species were living (Epstein et al., 1953). Surface seawater δ¹⁸O reflects the local precipitation-evaporation balance, local freshwater inputs like river runoff, and changes in global ice volume (Craig & Gordon, 1965; Shackleton & Opdyke, 1973). As for foraminifera Mg/Ca, biological and environmental factors influence foraminifera δ¹⁸O, including seasonality and depth habitat of the foraminifera, seawater pH, dissolution and post-depositional diagenesis (Pearson, 2012; Raymo et al., 2018).

In principle, the deep-sea δ¹⁸O<sub>benthic</sub> record is dependent on the δ¹⁸O composition of seawater and deep-sea temperature (Shackleton & Opdyke, 1973). However, in contrast to the surface (planktonic) signals, the δ¹⁸O<sub>benthic</sub> record is considered to have a strong influence of global ice volume and a more minor influence of local temperature and salinity changes, as deep waters are expected to be more uniform through time and space (Ahn et al., 2017; Lisiecki & Raymo, 2005; Rohling, 2013; Waelbroeck et al., 2002). When quantifying past ice volume change, the temperature component of the δ¹⁸O<sub>benthic</sub> record is assumed to have the largest uncertainty (e.g. Evans et al., 2016; Raymo et al., 2018). Here, we compare δ¹⁸O<sub>benthic</sub>
records from 22 sites which meet the PlioVAR criteria (Section 2.1), which span a range of regions and water depths to explore local/regional and global expressions of climate change. We did not include $\delta^{18}O_{benthic}$ records where more than one species had been used to generate the time series, so that we avoided any potential artificial introduction of shifts in the data when switching between species. Data are from *Cibicides* spp., *Cibicides wuellerstorfi*, *Uvigerina* spp., and *Cibicides mundulus* (Table 2).

To better understand the link between changes in ocean temperature and global ice volume, we examine indirect evidence for sea-level variability in the mPWP and early Pleistocene. Global sea-level variability is a measure of the magnitude and frequency of past ice volume fluctuations; direct records (e.g. paleo shorelines, speleothems, shallow marine sequences) are, however, spatially and temporally sparse (e.g. Miller et al., 2012 and references therein; Grant et al., 2019; Rovere et al., 2020). Age uncertainties associated with discontinuous records also makes it challenging to attribute sea-level maxima to single interglacials (e.g. Dumitru et al., 2019). To systematically assess sea-level variability through the mPWP and iNHG, we draw on the continuous but indirect sea-level estimates based on geochemical proxies from marine sediments (Miller et al., 2020; Rohling et al., 2021) and inverse modeling of $\delta^{18}O_{benthic}$ (Berends et al., 2021). Uncertainties come from the relative influence of temperature and diagenesis on $\delta^{18}O$ noted above, as well as the unconstrained $\delta^{18}O$ composition of polar ice sheets (Gasson et al., 2014) and varying ocean water masses. Here, we analyze three sea-level data sets generated using $\delta^{18}O_{benthic}$, which have been calibrated to sea level by accounting for the influence of bottom water temperature variability through time (Miller et al. (2020), hereafter “Miller2020”), by accounting for non-linear changes to ice sheet $\delta^{18}O$ during ice-sheet growth (Rohling et al. (2014), hereafter “Rohling2014”), or by undertaking an inverse modeling of temperature and ice volume contributions to $\delta^{18}O$ (Berends et al., 2021, hereafter “Berends2021”).

## 2. Materials and Methods

To address the questions posed in Section 1, the PAGES-PlioVAR group considered the window of interest for this study to span from 3.3 to 2.4 Ma. This definition was used to ensure our target window included the M2 glacial stage and the mPWP which immediately followed it, as well as the onset and intensification of glacial-interglacial cycles in $\delta^{18}O_{benthic}$ records which align with iNHG (Figure 1) (McClymont et al., 2017).

The majority of the sites we examined were recovered by the Deep Sea Drilling Project (DSDP), the Ocean Drilling Program (ODP), the Integrated Ocean Drilling Program and the International Ocean Discovery Program (both IODP) (Figure 2, Table 1). The specific location at which each sediment sequence was recovered is termed the Site, and it is IODP nomenclature to name the sequence by a combination of the drilling program and the site number. We use this nomenclature in our text, tables and graphics (e.g. ODP Site 907 refers to Site 907 recovered by the ODP).

### 2.1. Age Models

Marine sediment sequences can be placed on an age scale (a “chronology”) using a variety of methods. These include identifying events of known age, such as reversals in Earth's magnetic field and the first-appearance or extinction events in the fossil record (e.g., Gradstein et al., 2012). However, the intervals between these events can be hundreds of thousands of years apart. To explore climate changes on orbital timescales (see Section 1), we take advantage of the global impact of ice-sheet growth on the stable oxygen isotope ratio of seawater. The preferential storage of the $^{16}$O isotope in continental ice during glacial stages leads to an increase in seawater $\delta^{18}O$, which is then recorded in $\delta^{18}O_{benthic}$ (see also the Technical Box). Although individual $\delta^{18}O_{benthic}$ records will have some local influences, the approach of “stacking” (or averaging) across multiple sites has been shown to increase the signal-to-noise ratio of $\delta^{18}O_{benthic}$, and in turn provide a reference framework against which other sites can be aligned and then dated (Ahn et al., 2017; Lisiecki & Raymo, 2005). Regional differences in the expression of $\delta^{18}O_{benthic}$ and its relationship to orbital forcing have been determined (e.g. Caballero-Gill et al., 2019; Tian et al., 2002; Wilkens et al., 2017) which may hamper the accuracy of $\delta^{18}O_{benthic}$-derived stratigraphy. However, for
<caption>**Table 1**</caption>

The PlioVAR-Synthesized Records of Climate Change Spanning the Mid-Piacenzian Warm Period and Intensification of Northern Hemisphere Glaciation, Analyzed to Generate Figures 3–7 (Available at McClymont et al., 2023)

<table><thead><tr><th>Name</th><th>Lat.<br>(°N)</th><th>Long.<br>(°E)</th><th>Water<br>depth<br>(m)</th><th>Basin</th><th>Age model</th><th>Proxy</th><th>Original publication</th><th>Database source (URL or DOI)</th></tr></thead><tbody><tr><td>ODP 907</td><td>69.24</td><td>-12.7</td><td>1.802</td><td>natl</td><td>PlioVAR</td><td>U<sup>K</sup><sub>-37</sub></td><td>Clotten et al. (2018)</td><td>https://doi.org/10.1594/PANGAEA.877308</td></tr><tr><td>ODP 982</td><td>57.52</td><td>-15.87</td><td>1,134</td><td>natl</td><td>PlioVAR</td><td>U<sup>K</sup><sub>-37</sub></td><td>Lawrence et al. (2009)</td><td>https://doi.org/10.25921/52j8-r05</td></tr><tr><td>ODP 982</td><td>57.52</td><td>-15.87</td><td>1,134</td><td>natl</td><td>Lisiecki and<br>Raymo (2005)</td><td>δ<sup>18</sup>O benthic</td><td>Lisiecki and Raymo (2005)</td><td>https://doi.org/10.1594/PANGAEA.913056</td></tr><tr><td>UI313</td><td>41</td><td>-32.96</td><td>3,413</td><td>natl</td><td>Naafs et al. (2020)</td><td>U<sup>K</sup><sub>-37</sub></td><td>Naafs et al. (2010)</td><td>https://doi.org/10.1594/PANGAEA.913056</td></tr><tr><td>UI313</td><td>41</td><td>-32.96</td><td>3,413</td><td>natl</td><td>Naafs et al. (2020)</td><td>Mg/Ca (G.<br>bulloides)</td><td>Hennissen et al. (2014), De<br>Schepper et al. (2013)</td><td>https://doi.org/10.1594/PANGAEA.865414, https://doi.org/10.1594/PANGAEA.804675</td></tr><tr><td>UI313</td><td>41.00</td><td>-32.96</td><td>3,413</td><td>natl</td><td>Naafs et al. (2020)</td><td>δ<sup>18</sup>O benthic</td><td>Bolton et al. (2010)</td><td>https://doi.org/10.1594/PANGAEA.761444</td></tr><tr><td>UI313</td><td>41.00</td><td>-32.96</td><td>3,413</td><td>natl</td><td>Naafs et al. (2020)</td><td>δ<sup>18</sup>O plank.<br>(G.<br>bulloides)</td><td>Hennissen et al. (2014), De<br>Schepper et al. (2013)</td><td>https://doi.org/10.1594/PANGAEA.865414, https://doi.org/10.1594/PANGAEA.804675</td></tr><tr><td>ODP 999</td><td>12.74</td><td>-78.74</td><td>2,828</td><td>natl</td><td>Steph et al. (2010),<br>De Schepper<br>et al. (2013)</td><td>Mg/Ca (T.<br>sacculifer)</td><td>Groeneveld et al. (2014), De<br>Schepper et al. (2013)</td><td>https://doi.org/10.1594/PANGAEA.834675, https://doi.org/10.1594/PANGAEA.315652, https://doi.org/10.1594/PANGAEA.804671</td></tr><tr><td>ODP 999</td><td>12.74</td><td>-78.74</td><td>2,828</td><td>natl</td><td>De Schepper<br>et al. (2013)</td><td>δ<sup>18</sup>O benthic</td><td>De Schepper et al. (2013),<br>Haug and<br>Tiedemann (1998), and<br>Steph et al. (2010)</td><td>https://doi.org/10.1594/PANGAEA.804671, https://doi.org/10.1594/PANGAEA.789866, https://doi.org/10.1029/2008PA001645</td></tr><tr><td>ODP 999</td><td>12.74</td><td>-78.74</td><td>2,828</td><td>natl</td><td>De Schepper<br>et al. (2013)</td><td>δ<sup>18</sup>O plank.<br>(G.<br>sacculifer)</td><td>De Schepper et al. (2013),<br>Haug and<br>Tiedemann (1998)</td><td>https://doi.org/10.1594/PANGAEA.804671, https://doi.org/10.1594/PANGAEA.789867</td></tr><tr><td>ODP 662</td><td>-1.39</td><td>-11.74</td><td>3,814</td><td>natl</td><td>PlioVAR</td><td>U<sup>K</sup><sub>-37</sub></td><td>Herbert et al. (2010)</td><td>https://doi.org/10.1594/PANGAEA.874748</td></tr><tr><td>ODP 662</td><td>-1.39</td><td>-11.74</td><td>3,814</td><td>natl</td><td>Lisiecki and<br>Raymo (2005)</td><td>δ<sup>18</sup>O benthic</td><td>Lisiecki and Raymo (2005)</td><td>http://lorraine-lisiecki.com/stack.html</td></tr><tr><td>ODP 978</td><td>36.23</td><td>-2.06</td><td>1,930</td><td>natl</td><td>PlioVAR</td><td>U<sup>K</sup><sub>-37</sub></td><td>Khélifi et al. (2014)</td><td>https://doi.org/10.1594/PANGAEA.863960</td></tr><tr><td>ODP 978</td><td>36.23</td><td>-2.06</td><td>1,930</td><td>natl</td><td>PlioVAR</td><td>δ<sup>18</sup>O benthic</td><td>Khélifi et al. (2014)</td><td>https://doi.org/10.1594/PANGAEA.863938</td></tr><tr><td>ODP 978</td><td>36.23</td><td>-2.06</td><td>1,930</td><td>natl</td><td>PlioVAR</td><td>δ<sup>18</sup>O plank.<br>(G.<br>ruber)</td><td>García-Gallardo et al. (2018)</td><td>https://doi.org/10.5194/cp-14-339-2018</td></tr><tr><td>Sicily Punta<br>Piccola<br>("SPP")</td><td>37.29</td><td>13.49</td><td>NA</td><td>natl</td><td>Herbert et al. (2015)</td><td>U<sup>K</sup><sub>-37</sub></td><td>Herbert et al. (2015)</td><td>https://doi.org/10.25921/52j8-r05</td></tr><tr><td>ODP 704</td><td>-46.88</td><td>7.42</td><td>2.542</td><td>satl</td><td>Hodell and<br>Venz (1992)</td><td>δ<sup>18</sup>O plank.<br>(G.<br>bulloides)</td><td>Hodell and Venz (1992)</td><td>https://doi.org/10.1594/PANGAEA.58717</td></tr><tr><td>ODP 1082</td><td>-21.09</td><td>11.82</td><td>853</td><td>satl</td><td>Etourneau<br>et al. (2009)</td><td>U<sup>K</sup><sub>-37</sub></td><td>Etourneau et al. (2009)</td><td>https://doi.org/10.1594/PANGAEA.729403</td></tr><tr><td>ODP 1082</td><td>-21.09</td><td>11.82</td><td>853</td><td>satl</td><td>PlioVAR</td><td>δ<sup>18</sup>O benthic</td><td>Dupont et al. (2005)</td><td>https://doi.org/10.1594/PANGAEA.316972</td></tr></tbody></table>
<table><caption>**Table 1**
*Continued*</caption><thead><tr><th>Name</th><th>Lat.<br>(°N)</th><th>Long.<br>(°E)</th><th>Water<br>depth<br>(m)</th><th>Basin</th><th>Age model</th><th>Proxy</th><th>Original publication</th><th>Database source (URL or DOI)</th></tr></thead><tbody><tr><td>ODP 1087</td><td>-31.47</td><td>15.31</td><td>1,372</td><td>satl</td><td>Petrick et al. (2018)</td><td>U<sup>K</sup><sub>37</sub></td><td>Petrick et al. (2018)</td><td>https://doi.org/10.1594/PANGAEA.888605</td></tr><tr><td>ODP 1264</td><td>-28.53</td><td>2.85</td><td>2,507</td><td>satl</td><td>Bell et al. (2014)</td><td>δ<sup>18</sup>O benthic</td><td>Bell et al. (2014)</td><td>https://doi.org/10.25921/bk2a-2572</td></tr><tr><td>ODP 1267</td><td>-28.10</td><td>1.71</td><td>4,355</td><td>satl</td><td>Bell et al. (2014)</td><td>δ<sup>18</sup>O benthic</td><td>Bell et al. (2014)</td><td>https://doi.org/10.25921/bk2a-2572</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>113.29</td><td>2,772</td><td>pac</td><td>L. Li et al. (2011)</td><td>U<sup>K</sup><sub>37</sub></td><td>L. Li et al. (2011)</td><td>https://doi.org/10.1594/PANGAEA.956158</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>113.29</td><td>2,772</td><td>pac</td><td>Tian et al. (2006)</td><td>Mg/Ca (G.<br>ruber)</td><td>Tian et al. (2006)</td><td>https://doi.org/10.1594/PANGAEA.707839</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>113.29</td><td>2,772</td><td>pac</td><td>Tian et al. (2006)</td><td>δ<sup>18</sup>O benthic</td><td>Tian et al. (2006)</td><td>https://doi.org/10.1594/PANGAEA.700904</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>113.29</td><td>2,772</td><td>pac</td><td>Cheng et al. (2004)</td><td>δ<sup>18</sup>O plank.<br>(G.<br>ruber)</td><td>Cheng et al. (2004)</td><td>https://doi.org/10.1594/PANGAEA.784148</td></tr><tr><td>ODP 1148</td><td>18.84</td><td>116.57</td><td>3,294</td><td>pac</td><td>Jian et al. (2003)</td><td>δ<sup>18</sup>O benthic</td><td>Cheng et al. (2004)</td><td>https://doi.org/10.1594/PANGAEA.784180</td></tr><tr><td>ODP 1148</td><td>18.84</td><td>116.57</td><td>3,294</td><td>pac</td><td>Jian et al. (2003)</td><td>δ<sup>18</sup>O plank.<br>(G.<br>ruber)</td><td>Cheng et al. (2004)</td><td>https://doi.org/10.1594/PANGAEA.784180</td></tr><tr><td>ODP 1208</td><td>36.13</td><td>158.20</td><td>3,346</td><td>pac</td><td>Venti and<br>Billups (2012)</td><td>U<sup>K</sup><sub>37</sub></td><td>Abell et al. (2021), Venti<br>et al. (2013)</td><td>https://figshare.com/articles/dataset/Dust_SST_and_Productivity_Data_<br>for_ODP_1208A_and_ODP_885_886/12472646</td></tr><tr><td>ODP 1208</td><td>36.13</td><td>158.20</td><td>3,346</td><td>pac</td><td>Venti and<br>Billups (2012)</td><td>δ<sup>18</sup>O benthic</td><td>Venti and Billups (2012)</td><td>https://doi.org/10.25921/hamy-bb98</td></tr><tr><td>ODP 1241</td><td>5.84</td><td>-86.44</td><td>2,027</td><td>pac</td><td>Tiedemann<br>et al. (2007)</td><td>Mg/Ca (T.<br>sacculifer)</td><td>Groeneveld and<br>Tiedemann (2005)</td><td>https://doi.org/10.1594/PANGAEA.315654</td></tr><tr><td>ODP 1241</td><td>5.84</td><td>-86.44</td><td>2,027</td><td>pac</td><td>Tiedemann<br>et al. (2007)</td><td>δ<sup>18</sup>O benthic</td><td>Tiedemann et al. (2007)</td><td>https://doi.org/10.1594/PANGAEA.774009</td></tr><tr><td>ODP 1012</td><td>32.28</td><td>-118.38</td><td>1,772</td><td>pac</td><td>Brierley et al. (2009)</td><td>U<sup>K</sup><sub>37</sub></td><td>Brierley et al. (2009)</td><td>https://doi.org/10.1594/PANGAEA.956158</td></tr><tr><td>U1417</td><td>56.96</td><td>-147.11</td><td>4,187</td><td>pac</td><td>Sánchez-Montes<br>et al. (2020)</td><td>U<sup>K</sup><sub>37</sub></td><td>Sánchez-Montes et al. (2020)</td><td>https://doi.org/10.1594/PANGAEA.899064</td></tr><tr><td>ODP 806</td><td>0.32</td><td>159.36</td><td>2,520</td><td>pac</td><td>PlioVAR</td><td>Mg/Ca (T.<br>sacculifer)</td><td>Wara et al. (2005)</td><td>https://www.science.org/action/<br>downloadSupplement?doi=10.1126/%2Fscience.1112596&file=wara.<br>som.rev1.pdf</td></tr><tr><td>ODP 806</td><td>0.32</td><td>159.36</td><td>2,520</td><td>pac</td><td>PlioVAR</td><td>δ<sup>18</sup>O benthic</td><td>Karas et al. (2009)</td><td>https://doi.org/10.1038/ngeo0520</td></tr><tr><td>ODP 846</td><td>-3.09</td><td>-90.82</td><td>3,296</td><td>pac</td><td>Lawrence<br>et al. (2006)</td><td>U<sup>K</sup><sub>37</sub></td><td>Lawrence et al. (2006)</td><td>https://doi.org/10.25921/54bt-b232</td></tr><tr><td>ODP 846</td><td>-3.09</td><td>-90.82</td><td>3,296</td><td>pac</td><td>Shackleton<br>et al. (1995)</td><td>δ<sup>18</sup>O benthic</td><td>Shackleton et al. (1995)</td><td>https://doi.org/10.1594/PANGAEA.696450</td></tr><tr><td>ODP 849</td><td>0.18</td><td>-110.52</td><td>3,839</td><td>pac</td><td>Mix et al. (1995)</td><td>δ<sup>18</sup>O benthic</td><td>Mix et al. (1995)</td><td>https://doi.org/10.1594/PANGAEA.701400</td></tr><tr><td>ODP 214</td><td>-11.30</td><td>88.70</td><td>1,665</td><td>ind</td><td>Karas et al. (2009)</td><td>Mg/Ca (G.<br>ruber)</td><td>Karas et al. (2009)</td><td>https://doi.org/10.25921/b227-3f44</td></tr><tr><td>ODP 214</td><td>-11.30</td><td>88.70</td><td>1,665</td><td>ind</td><td>Karas et al. (2009)</td><td>δ<sup>18</sup>O benthic</td><td>Karas et al. (2009)</td><td>https://doi.org/10.25921/b227-3f44</td></tr></tbody></table>
<table><caption>Table 1<br>Continued</caption><thead><tr><th>Name</th><th>Lat.<br>(°N)</th><th>Long.<br>(°E)</th><th>Water<br>depth<br>(m)</th><th>Basin</th><th>Age model</th><th>Proxy</th><th>Original publication</th><th>Database source (URL or DOI)</th></tr></thead><tbody><tr><td>ODP 214</td><td>-11.30</td><td>88.70</td><td>1.665</td><td>ind</td><td>Karas et al. (2009)</td><td>δ<sup>18</sup>O plank.<br>(G. ruber)</td><td>Karas et al. (2009)</td><td>https://doi.org/10.25921/b2z7-3f44</td></tr><tr><td>ODP 722</td><td>16.60</td><td>59.80</td><td>2.022</td><td>ind</td><td>Herbert et al. (2010)</td><td>U<sup>K</sup>'<sub>37</sub></td><td>Herbert et al. (2010)</td><td>https://doi.org/10.25921/we6n-qe21</td></tr><tr><td>ODP 722</td><td>16.60</td><td>59.80</td><td>2.022</td><td>ind</td><td>Clemens et al. (1996)</td><td>δ<sup>18</sup>O benthic</td><td>Clemens et al. (1996)</td><td>https://doi.org/10.25921/f8rt-w842</td></tr><tr><td>ODP 758</td><td>5.40</td><td>90.40</td><td>2.925</td><td>ind</td><td>Chen et al. (1995)</td><td>δ<sup>18</sup>O benthic</td><td>Chen et al. (1995)</td><td>https://doi.org/10.1594/PANGAEA.696412</td></tr><tr><td>U1463</td><td>-18.97</td><td>117.62</td><td>145</td><td>ind</td><td>Smith and<br>Castañeda (2020)</td><td>TEX<sub>96</sub></td><td>Smith and Castañeda (2020)</td><td>https://doi.org/10.25921/e1p0-9i22</td></tr><tr><td>DSDP 593</td><td>-40.51</td><td>167.67</td><td>1.088</td><td>so</td><td>McClymont<br>et al. (2016)</td><td>U<sup>K</sup>'<sub>37</sub></td><td>McClymont et al. (2016)</td><td>https://doi.org/10.25921/13rs-tp69</td></tr><tr><td>DSDP 593</td><td>-40.51</td><td>167.67</td><td>1.088</td><td>so</td><td>McClymont<br>et al. (2016)</td><td>δ<sup>18</sup>O benthic</td><td>McClymont et al. (2016)</td><td>https://doi.org/10.25921/13rs-tp69</td></tr><tr><td>DSDP 594</td><td>-45.68</td><td>174.96</td><td>1.204</td><td>so</td><td>Caballero-Gill<br>et al. (2019)</td><td>U<sup>K</sup>'<sub>37</sub></td><td>Caballero-Gill et al. (2019)</td><td>https://doi.org/10.1594/PANGAEA.898167</td></tr><tr><td>DSDP 594</td><td>-45.68</td><td>174.96</td><td>1.204</td><td>so</td><td>Caballero-Gill<br>et al. (2019)</td><td>δ<sup>18</sup>O benthic</td><td>Caballero-Gill et al. (2019)</td><td>https://doi.org/10.1594/PANGAEA.898167</td></tr><tr><td>ODP 1088</td><td>-41.14</td><td>13.56</td><td>2.082</td><td>so</td><td>Hodell and Venz-<br>Curtis (2006)</td><td>δ<sup>18</sup>O benthic</td><td>Hodell and<br>Venz-Curtis (2006)</td><td>https://doi.org/10.25921/c3ah-c727</td></tr><tr><td>ODP 1090</td><td>-42.91</td><td>8.90</td><td>3.699</td><td>so</td><td>Martinez-Garcia<br>et al. (2011)</td><td>U<sup>K</sup>'<sub>37</sub></td><td>Martinez-Garcia et al. (2010)</td><td>https://doi.org/10.1594/PANGAEA.771708</td></tr><tr><td>ODP 1092</td><td>-46.41</td><td>7.08</td><td>1.974</td><td>so</td><td>Andersson<br>et al. (2002)</td><td>δ<sup>18</sup>O benthic</td><td>Andersson et al. (2002)</td><td>https://doi.org/10.1594/PANGAEA.956158</td></tr></tbody><tfoot><tr><td colspan="10">Note. Where the age model is noted as “PlioVAR” the record has been updated after the original publication. Basin refers to North Atlantic (natl), South Atlantic (satl), Pacific (pac), Indian (ind) and Southern Oceans (so).</td></tr></tfoot></table>
<table><caption>**Table 2**

Sites Used for δ¹⁸O<sub>benthic</sub> Probability Distribution Analysis</caption><thead><tr><th>Site</th><th>Benthic species</th><th>Latitude (°N)</th><th>Longitude (°E)</th><th>Water depth (m)</th><th>Ocean basin</th><th>Period</th></tr></thead><tbody><tr><td>U1313</td><td>C. wuellerstorfi</td><td>41.00</td><td>-32.96</td><td>3,426</td><td>North Atlantic</td><td>3.3-2.4</td></tr><tr><td>ODP 978*</td><td>C. wuellerstorfi</td><td>36.23</td><td>-2.06</td><td>1,940</td><td>Eastern Mediterranean</td><td>3.3-2.7</td></tr><tr><td>ODP 1208</td><td>C. wuellerstorfi</td><td>36.13</td><td>158.20</td><td>3,346</td><td>Northwest Pacific</td><td>3.3-2.4</td></tr><tr><td>ODP 999</td><td>C. wuellerstorfi</td><td>12.74</td><td>-78.74</td><td>2,838</td><td>Mid Atlantic (Caribbean Sea)</td><td>3.3-2.4</td></tr><tr><td>ODP 1143</td><td>C. wuellerstorfi</td><td>9.36</td><td>113.29</td><td>2,772</td><td>South China Sea</td><td>3.3-2.4</td></tr><tr><td>ODP 1241</td><td>C. wuellerstorfi</td><td>5.84</td><td>-86.44</td><td>2,027</td><td>Eastern Equatorial Pacific</td><td>3.3-2.4</td></tr><tr><td>ODP 758</td><td>C. wuellerstorfi</td><td>5.38</td><td>90.36</td><td>2,923</td><td>Mid Indian Ocean</td><td>3.3-2.4</td></tr><tr><td>ODP 806</td><td>C. wuellerstorfi</td><td>0.32</td><td>159.36</td><td>2,520</td><td>Western Equatorial Pacific</td><td>3.3-2.4</td></tr><tr><td>ODP 849</td><td>C. wuellerstorfi</td><td>0.18</td><td>-110.52</td><td>3,850</td><td>Eastern Equatorial Pacific</td><td>3.3-2.4</td></tr><tr><td>ODP 662</td><td>Cibs. spp.</td><td>-1.39</td><td>-11.74</td><td>3,821</td><td>Mid Atlantic Ocean</td><td>3.3-2.4</td></tr><tr><td>ODP 849</td><td>C. wuellerstorfi</td><td>0.18</td><td>-110.52</td><td>3,839</td><td>Eastern Equatorial Pacific</td><td>3.3-2.4</td></tr><tr><td>ODP 1082</td><td>C. wuellerstorfi</td><td>-21.09</td><td>11.82</td><td>1,290</td><td>South Atlantic</td><td>3.3-2.4</td></tr><tr><td>ODP 1267</td><td>C. wuellerstorfi</td><td>-28.10</td><td>1.71</td><td>4,356</td><td>South Atlantic</td><td>2.7-2.4</td></tr><tr><td>ODP 1264</td><td>C. wuellerstorfi</td><td>-28.53</td><td>2.85</td><td>2,507</td><td>South Atlantic</td><td>2.7-2.4</td></tr><tr><td>ODP 1088</td><td>Cibs. spp.</td><td>-41.14</td><td>13.56</td><td>2,254</td><td>South Atlantic</td><td>3.3-2.4</td></tr><tr><td>DSDP 594*</td><td>Cibs. spp.</td><td>-45.52</td><td>174.95</td><td>1,204</td><td>West Pacific</td><td>3.3-2.4</td></tr></tbody></table>

Note. Sites are ordered by latitude as shown in Figure 6, and mapped on Figure 2. Sites which do not have data spanning the whole interval of interest (3.3–2.4 Ma) are asterisked (ODP Site 978, DSDP Site 594).

our interval of study, independent high-resolution δ¹⁸O<sub>benthic</sub> stratigraphies from the equatorial Atlantic confirm the LR04 age assignments (Wilkens et al., 2017).

Following protocols established by PAGES-PlioVAR (McClymont et al., 2017, 2020), sites were only included in this synthesis if they had: (a) δ¹⁸O<sub>benthic</sub> data at ≤10 kyr resolution, tied to the global benthic stacks (LR04: Lisiecki & Raymo, 2005; Prob-stack: Ahn et al., 2017), and/or (b) climate data at ≤10 kyr resolution, constrained by paleomagnetic tie-points including the Gauss/Matuyama boundary (2.581 Ma), upper Mammoth (3.207 Ma) and lower Mammoth (3.330 Ma) (Gradstein et al., 2012). The age controls for all 32 sites which met these thresholds were reviewed, and the age models for several sites were revised compared to their published chronologies (Table 1 and Figure S1 in Supporting Information S1). Two sites were nevertheless excluded because although they met the age model requirements outlined here, the δ¹⁸O<sub>benthic</sub> data were unable to discern glacial-interglacial cycling as identified in the global stacks (n = 2; Figure S1 in Supporting Information S1). Fifty-three climate records met our criteria outlined above (the “PlioVAR data sets”), which include 5 climate proxies across 32 sites (Figure 1, Table 1).

## 2.2. Marine Proxies for Ocean Temperature, Ice Volume and Sea Level

A wide range of principles underpins the methods we employ here to reconstruct climate variables. The technical details and the considerations we made in our selections are provided in the accompanying Technical Box. In brief, we synthesized SST reconstructions generated by three proxies (UK' 37, TEX<sub>86</sub>, and Mg/Ca). The UK' 37 index is a ratio of organic molecules (alkenones) synthesized by selected haptophyte algae, which photosynthesize in the surface ocean (Marlowe et al., 1984; Prahl & Wakeham, 1987). The TEX<sub>86</sub> index is a ratio of organic molecules (isoprenoidal glycerol dialkyl glycerol tetraethers, isoGDGTs) synthesized by selected marine Archaea, in particular the ammonia oxidizing *Thaumarchaeota* (Schouten et al., 2002, 2013). SSTs are also reconstructed from the Mg/Ca ratio in the shells of planktonic foraminifera: single-celled protists found at a range of depths in the upper parts of the ocean water column. As each of these three SST proxies has different biological sources and environmental influences, there is the potential to extract a rich suite of information by drawing on the similarities and differences of the reconstructed SSTs (discussed further in Section 3).
We also synthesized foraminiferal $\delta^{18}$O records from species living in the upper water column (i.e., near the surface, $\delta^{18}$O<sub>planktonic</sub>) and at the sea floor ($\delta^{18}$O<sub>benthic</sub>). Foraminifera $\delta^{18}$O record both the influences of the temperature and $\delta^{18}$O composition of the seawater in which they live (Epstein et al., 1953). The seawater $\delta^{18}$O is also influenced by temperature and salinity, especially in the surface ocean, whereas in the deep sea an increasing influence of global ice volume is recognized (Shackleton & Opdyke, 1973) but not always easy to quantify (Evans et al., 2016; Raymo et al., 2018). The properties of the intermediate and deep waters which bathe the sites thus reflect surface-ocean conditions in the regions where they were formed through density-driven convection; in our compilation these include the high latitude oceans of both hemispheres and the Mediterranean Sea (e.g. Hodell & Venz-Curtis, 2006; Khélifi et al., 2014). As a result, multi-proxy analyses from a single site have the potential to record asynchronous changes between the surface and deep ocean. We also analyzed the results of several studies which reconstructed the sea-level variability contributions to $\delta^{18}$O<sub>benthic</sub> by accounting for differences in bottom water temperature variability through time (Miller et al., 2020), by accounting for non-linear changes to ice sheet $\delta^{18}$O during ice-sheet growth (Rohling et al., 2014), and by undertaking an inverse modeling of temperature and ice volume contributions to $\delta^{18}$O<sub>benthic</sub> (Berends et al., 2021).

## 2.3. Statistical Analysis

We adopted several statistical approaches to assess the character of the mPWP and to investigate how the climate system evolved across the iNHG. Previous analysis of the 5–0 Ma LR04 $\delta^{18}$O<sub>benthic</sub> stack using a Bayesian Change Point algorithm identified a change point at 2.73 Ma (±0.1 Ma), closely aligning with iNHG, which was marked by an increase in the obliquity signal (Ruggieri, 2013). To test whether regime boundaries also existed in the PlioVAR data sets (e.g. a shift in the mean or variance of the climate data), we applied the same Bayesian Change Point algorithm to the 53 individual proxy datasets and the three sea-level records which had high-resolution data for the whole 3.3–2.4 Ma interval (Table 1). This algorithm calculates the probability density of the input time series, and uses a regression model and the combined application of forward recursion and Bayes rule to compute change point locations and their probability. As input, the program requires the number of change points to be found, the minimum distance between change points, and three hyperparameters that characterize the prior distribution of regression coefficients and residual variance (Ruggieri, 2013). The hyperparameters are used to scale the prior distributions of regression coefficients (based on a multivariate normal distribution) and residual variance (scaled-inverse $\chi^2$ distribution). Analyses using this program were performed using the same hyperparameters as used by Ruggieri (2013) for the Plio-Pleistocene LR04 $\delta^{18}$O<sub>benthic</sub> stack.

As we sought to identify whether iNHG represented a single shift in climate regime, as shown by Ruggieri (2013), we set the input to find only one viable change point within the analysis interval. This change point represents the mid-point of the climate shift that is, the point at which the preceding climate is different to the one that follows, rather than identifying the onset of a potential transition (Figure S3 in Supporting Information S1). To avoid single large glacial or interglacial peaks or data outliers artificially being indicated as change points, and accounting for the temporal resolution of our data sets, we fixed the minimum distance of detection to be equivalent to ~10,000 years, regardless of the sedimentation rates at each site. Furthermore, some discontinuous records (8 out of 56) were interpolated using a cubic spline to their median time step prior to analysis.

We also assessed and compared the amplitude of glacial-interglacial variability spanning three periods for the $\delta^{18}$O<sub>benthic</sub> and SST data. We refer to these as intervals P-1, P-2, and P-3 throughout the manuscript (Figure 1). Interval P-1 broadly aligns with the mPWP (3.3–3.0 Ma). Interval P-2 (3.0–2.7 Ma) represents the transition toward Interval P-3 (2.7–2.4 Ma), which is marked by the onset of consistently elevated IRD deposition in the high northern latitudes (Figure 1, Section 1). By comparing the amplitudes of glacial-interglacial variability we reduce the reliance on peak interglacial correlation and influence of individual age model uncertainty which might have impacted our comparison of globally distributed records with varying temporal resolution. For each record in the three intervals (P-1, P-2, P-3), we also generated the probability distribution of glacial-interglacial amplitudes for selected $\delta^{18}$O<sub>benthic</sub> and SST records (after Grant & Naish, 2021). To capture glacial-interglacial variability, the ranges of amplitude (maximum-minimum values) were calculated within a centered 41 kyr moving window, selected as the primary periodicity modulating glacial cyclicity, using
R package Astrochron (Meyers et al., 2021). A minimum sampling of <10 kyr is required to resolve peak glacial-interglacial amplitudes. Data sets which met this requirement (n = 49) were linearly interpolated to the mean sampling resolution for each site (1–5 kyr), which was also used for the window time-step (i.e., the time between centered windows).

Finally, we assess the amplitude of glacial-interglacial sea-level cycles for the same three time intervals following Grant and Naish (2021). The full amplitude of sea-level change from colder-than-present glacial states to warmer interglacials provides insight on ice-sheet sensitivity, the magnitude of glacial-interglacial variability, and also allowed us to evaluate how changes in the cryosphere during the iNHG compared to the evolution of SSTs discussed in this manuscript. As we did for the δ¹⁸O_benthic and SST records, the glacial-interglacial amplitudes were generated as the maximum range within a 41 kyr moving-window, for intervals P-1, P-2, and P-3, using a mean sampling step of the record (2 kyr). The amplitudes derived from this approach are on a floating scale (i.e., they sit unanchored to a Holocene reference).

### 3. Are There Proxy-Specific Differences in the Records of SST Change Through Time? (Q1)

Previous global-scale syntheses of Pliocene climate have used reconstructed SSTs as key variables, especially for data-model comparison (e.g., Dowsett et al., 2012; Haywood et al., 2020; McClymont et al., 2020; Z. Zhang et al., 2021). This previous work has allowed quantification of climate response to pCO₂ forcing (Martinez-Boti et al., 2015), and enabled comparison to other geological time intervals where forcings and feedbacks may have been different (e.g. Pre-Industrial, Last Glacial Maximum, Eocene) (Fischer et al., 2018; Gulev et al., 2021). Although there is often good correspondence between different SST proxies, differences in absolute values and trends can arise which may vary by location but also through time (Lawrence & Woodard, 2017; McClymont et al., 2020). In turn, it can become challenging to calculate site-specific and global-scale temperature anomalies, and to undertake data-model comparison (e.g. Inglis et al., 2020; McClymont et al., 2020). Since each proxy system has a different biological origin with different biological and environmental controls on how the temperature signal is recorded (Section 2), we propose that there is the potential to understand the evolution of Pliocene and early Pleistocene environments more thoroughly by adopting a multi-proxy approach.

In this section, we examine SST records from three regional settings for which the most oNHG and iNHG temperature data have been generated: (a) the North Atlantic Ocean; (b) low latitude warm pools; (c) low latitude upwelling regions (Figure 2). The different oceanographic settings allow us to explore likely controls over the proxy signals. Since we have already accounted for a modest Mg/Ca<sub>seawater</sub> correction (Technical Box), offsets between foraminifera Mg/Ca reconstructions and other proxies need to be explained by other factors, which could include for example, the preferred season of growth, or the water depth where the source organisms were living (Bova et al., 2021). We sought to identify whether overlaps or offsets existed between data sets which might be consistent with our understanding of regional processes and biological controls. This was achieved by: (a) examining multi-proxy data from individual sites, where more than one proxy has been analyzed from the same sediment sequence; (b) examining SST data at a regional scale (Figure 3); and (c) comparing the PlioVAR data sets to published lower-resolution data. To facilitate comparison between sites where there may be different temporal resolutions of data, we plot PlioVAR SSTs as 20 kyr means centered on each glacial and interglacial maximum/minimum as determined by δ¹⁸O_benthic (Figure 3). SST maxima or minima do not always align with the peak glacial or interglacial states as determined by the marine isotope stages (MIS) of the δ¹⁸O_benthic record (e.g. Capron et al., 2014). However, our approach here is to use δ¹⁸O_benthic as a simple framework for identifying mean SSTs within broad ~20 kyr windows, to assess whether any proxy-specific bias occurs across iNHG.

#### 3.1. A Comparison of Multi-Proxy Temperature Records From Single Sites

An important limitation to our analysis is that only 2 of the 22 PlioVAR sites examined here have both high temporal resolution and multi-proxy SST data recovered from the same location to enable direct proxy-proxy comparison across the whole interval of study (3.3–2.4 Ma): at the North Atlantic site IODP U1313, and the South China Sea site ODP 1143 (Table 1; Figure 2).
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Regional comparisons of glacial and interglacial means for selected marine isotope stages (MIS) spanning the late Pliocene and early Pleistocene. Only the MIS with the highest frequency of sea-surface temperature (SST) data and where the glacial or interglacial duration was ≥20 kyr were selected. Data are presented as 20 kyr means centered on the glacial or interglacial maximum as defined by Lisiecki and Raymo (2005), with the 2σ range of SSTs within that 20 kyr interval shown by the vertical error bars. For data sources see Table 1, for site locations see Figure 2. (a) δ¹⁸O_benthic stack (Lisiecki & Raymo, 2005) and the MIS analyzed in panels (a–c) where even numbers (blue shading) highlight glacial stages, and odd numbers (orange shading) highlight interglacial (warm) stages. (b) North Atlantic sites. Note that the Site SPP (Mediterranean Sea) lies outside the North Atlantic but its SSTs are likely to be closely connected to changes there; (c) Low-latitude sites which are unaffected by upwelling (“warm pools”); (d) Low-latitude upwelling sites; (e) Temperatures are reconstructed using U<sub>37</sub><sup>K'</sup> (10 sites, squares) and foraminifera Mg/Ca (5 sites, triangles).</figcaption></figure>

Mg/Ca in *G. bulloides* consistently yields 1–2°C higher SSTs than U<sub>37</sub><sup>K'</sup> at the North Atlantic site (IODP Site U1313), excluding MIS 100 (~2.52 Ma) where Mg/Ca *G. bulloides* temperatures are within error of, but slightly cooler than, U<sub>37</sub><sup>K'</sup> (Figure 3a). At this latitude (41°N), U<sub>37</sub><sup>K'</sup> is calibrated to mean annual SST (Müller et al., 1998; Tierney & Tingley, 2018), suggesting that the *G. bulloides* temperatures may be reflecting warmer seasonal SSTs. We did not analyze the *G. ruber* Mg/Ca data from IODP Site U1313 (Bolton et al., 2010, 2018) due to discontinuous sampling across our study interval, but the reconstructed *G. ruber* SSTs are warmer than from *G*.
bulloides, consistent with modern summer production in *G. ruber* and spring production in *G. bulloides*, as also observed in other North Atlantic records (Friedrich et al., 2013; Hennissen et al., 2014; Robinson et al., 2008). The enhanced cooling in *G. bulloides* Mg/Ca during glacial stages (Figure 3a) thus indicates either enhanced spring cooling associated with iNHG, or a shift to deep (colder) growth temperatures, given that *G. bulloides* occupies a mixed-layer habitat (the upper ~60 m) (e.g. Schiebel et al., 1997) compared to the upper 10 m for $U_{37}^{K'}$ (Müller et al., 1998).

In contrast, in the low latitude warm pools (Figure 2), $U_{37}^{K'}$ SSTs are within proxy calibration uncertainty of the Mg/Ca temperature estimates, as shown in the South China Sea (ODP Site 1143) with Mg/Ca in *G. ruber* (Figure 3), but also when comparing a lower time resolution $U_{37}^{K'}$ record to the high-resolution *T. sacculifer* Mg/Ca estimates from the Caribbean Sea (ODP Site 999, Badger et al., 2013; Figure S4 in Supporting Information S1). However, $U_{37}^{K'}$ temperatures are consistently higher than Mg/Ca in the South China Sea. Strong seasonal variations in SST and mixed-layer depth are recorded today in the South China Sea, linked to the East Asian monsoon system (Twigt et al., 2007). A summer signal from *G. ruber* (Lin & Hsieh, 2007) is unable to explain the cooler Mg/Ca temperatures we reconstruct relative to the (mean annual) $U_{37}^{K'}$ signal; rather, the high winter fluxes of *G. ruber* observed today (Lin et al., 2011) may also have imparted a relatively cool signal during the late Pliocene and early Pleistocene. However, the reconstructed SSTs at ODP Site 1143 sit in the warmest part of the $U_{37}^{K'}$ calibration, where errors approach 4.4°C (Tierney & Tingley, 2018), meaning that in effect there is agreement, within error, between the two proxies in these high resolution records.

## 3.2. Regional-Scale Comparisons of Ocean Temperature Change

As an alternative to multi-proxy comparisons within individual sites, we compared the regional-scale glacial and interglacial patterns of SSTs from the PlioVAR data sets (Figure 3) to published lower-resolution data sets. In the North Atlantic Ocean, there is greater variability in $U_{37}^{K'}$ records between sites, than between the co-registered Mg/Ca and $U_{37}^{K'}$ data from IODP Site U1313 (Figure 3a). This variability demonstrates the sensitivity of these sites to changes in the position of the subtropical and subpolar gyres (e.g. IODP Site U1313, ODP Site 982) and the extent of Arctic surface water masses (e.g. ODP Site 907). With only one site containing both Mg/Ca and $U_{37}^{K'}$ data we are unable to assess whether there are any systematic differences in how SST is recorded by these proxies (e.g. season or depth of production). Several low-resolution multi-proxy time series are available from the North Atlantic Ocean, which did not meet the PlioVAR thresholds for analysis (Section 2.1; Figure 2), but they offer some insights into regional similarities and differences between our SST proxies.

For example, a reconstruction centered on MIS M2 (spanning MIS MG2-M1, ~3.34–3.24 Ma) showed comparable $U_{37}^{K'}$ trends among three sites along the path of the North Atlantic Current (spanning~40 to ~53°N) (De Schepper et al., 2013). SSTs from $U_{37}^{K'}$ were consistently ~2°C warmer than Mg/Ca-temperatures from *G. inflata* (De Schepper et al., 2013) consistent with the interpretation that *G. inflata* reflect deeper (~200–300 m) and thus cooler temperatures (De Schepper et al., 2013). Low-resolution $U_{37}^{K'}$ and $\text{TEX}_{86}$ temperatures ~53°N (DSDP Site 610) had no consistent offset through time, although their SSTs were within analytical and calibration errors (Naafs et al., 2020). In contrast, the $U_{37}^{K'}$ and $\text{TEX}_{86}$ temperatures were consistently warmer than *G. bulloides* Mg/Ca data from DSDP Site 610, where both proxies spanned the iNHG (Hennissen et al., 2014; Naafs et al., 2020). Cooler SSTs from *G. bulloides* might also reflect a spring or mixed-layer temperature signature here, as considered above for IODP Site U1313 (Friedrich et al., 2013; Hennissen et al., 2014; Robinson et al., 2008). Our regional-scale proxy comparison for the North Atlantic Ocean supports published interpretations of $U_{37}^{K'}$ and $\text{TEX}_{86}$ as proxies for mean annual or summer SSTs, compared to summer (*G. ruber*), spring or mixed layer (*G. bulloides*) or thermocline (*G. inflata*) proxies from foraminiferal Mg/Ca.

In the low latitude warm pools no long-term trends in glacial or interglacial SSTs were determined, regardless of the proxy (Figure 3b). However, both $U_{37}^{K'}$ and Mg/Ca SSTs in the South China Sea (ODP Site 1143) consistently exceed the Mg/Ca SSTs from all other warm pool sites, despite sitting more distal to the heart of the west Pacific warm pool (as recorded by proxy data from ODP Site 806, Figure 3b). The apparently warm Mg/Ca temperatures in the South China Sea may be explained by changes in salinity: a climate model which reproduces the lower zonal and meridional SST gradients shown in Pliocene proxy data (Burls et al., 2017) also generates saltier conditions at ODP Site 806 (+0.5 psu) and much saltier conditions at ODP Site 1143 (+2.5 psu). Given the salinity effect on
Mg/Ca (+5.4% per psu) (Gray & Evans, 2019), increased Pliocene salinity would cause the Mg/Ca temperatures in the South China Sea to appear ~1.3°C warmer than at ODP Site 806; without this effect, mPWP temperatures at the two sites are within ~0.3°C (Figure 3b). A difficulty with this explanation is that the low-resolution TEX86 temperatures from ODP Site 1143 are even warmer than both $U_{37}^{K'}$ and T. sacculifer Mg/Ca temperatures in the Pliocene (O'Brien et al., 2014), but cannot be explained by a response to elevated salinity alone. In contrast, Pliocene SSTs in the Caribbean Sea (ODP Site 999) are ~2.5°C cooler than the other warm pool sites (Figure 3b). Paleo-salinity changes are predicted to have been negligible (Burls et al., 2017; Rosenbloom et al., 2013), but seafloor dissolution may have been different (Groeneveld et al., 2014; Haug & Tiedemann, 1998); more corrosive deep waters in the Pliocene could have biased Mg/Ca-based temperatures to artificially cold values, which could be tested in the future at these sites by generating proxy data for bottom water carbonate ion saturation in tandem with Mg/Ca analysis (Rosenthal et al., 2022).

Although long-term evolution of the major upwelling systems have been intensively studied through the Pliocene and Pleistocene (e.g., Dekens et al., 2007; Etourneau et al., 2009; Rosell-Melé et al., 2014), the PlioVAR data sets are solely based on $U_{37}^{K'}$ analyses (Figure 3; Table 1). Yet, modern studies and multi-proxy analysis indicates the potential for seasonality or water depth to influence both the organic and inorganic proxies that we used (e.g. Hertzberg et al., 2016; Ho et al., 2011; Leduc et al., 2014; Lopes dos Santos et al., 2010; McClymont et al., 2012; Petrick et al., 2018; White & Ravelo, 2020). Additional low-resolution multi-proxy data are available at two upwelling sites from the south-east Atlantic, ODP Site 1087 (TEX86) and ODP Site 1082 (G. bulloides Mg/Ca), where $U_{37}^{K'}$ consistently records the warmest SSTs. The relatively cold TEX86 temperatures have been attributed to sub-surface production (Petrick et al., 2018; Seki et al., 2012), as is also observed in both late Pleistocene and modern upwelling systems (e.g. Hertzberg et al., 2016; Lopes dos Santos et al., 2010; McClymont et al., 2012) and discussed in detail for the Pliocene (Ford et al., 2015; White & Ravelo, 2020). Pliocene and Pleistocene SSTs from Mg/Ca in G. bulloides are consistently colder than $U_{37}^{K'}$ SSTs at ODP Site 1082, regardless of whether a Mg/Ca seawater correction is applied (up to 10°C, Leduc et al., 2014; ~1.5°C with Mg/Ca seawater correction, this study, Figure S2 in Supporting Information S1). The $U_{37}^{K'}$-Mg/Ca offset was attributed to G. bulloides recording a deeper, winter signal during upwelling intensification at ODP Site 1082 (Leduc et al., 2014). In turn, if $U_{37}^{K'}$ is recording warm-season SSTs (Leduc et al., 2014), this may partially explain some of the observed data-model offsets in the mPWP, given that the model results were of mean annual SST (McClymont et al., 2020). Care is thus needed at upwelling sites to select the most appropriate SST reconstructions (e.g whether annual or seasonal) and interpret their results in relation to both complex local oceanographic change and the imprint of global climate transitions.

### 3.3. A Comparison of Ocean Temperature Proxy Records: Summary and Recommendations

Our synthesis and proxy comparison of regional SST data sets spanning the mPWP and iNHG demonstrates that for most sites, $U_{37}^{K'}$, TEX86 and Mg/Ca in surface-dwelling foraminifera can provide robust reconstructions of past mean annual SST. We have also shown that under certain circumstances, local or regional influences over the season or depth of production by the source organisms can influence the temperature signals, which can be identified through multi-proxy analysis from individual sites.

Our analysis highlighted that care should be taken to ensure that G. bulloides Mg/Ca data is reflecting a surface-dwelling or annual signal when deciding whether to include it into global-scale syntheses of mean annual SSTs. This is because G. bulloides-derived SSTs may be biased toward spring or subsurface temperatures, but with a relationship that may not be constant through time, for example, if there have been orbital or longer-term changes in seasonality, mixed-layer depth, or upwelling location/intensity. TEX86 data may also show subsurface bias, and should be excluded from global SST syntheses where there is evidence for an increase in upwelling/ subsurface export, either directly (i.e., high GDGT-2/GDGT-3 ratios; van der Weijst et al., 2022) or indirectly (e.g. where microfossil assemblages highlight changing upwelling intensity). Incorporating seasonal or subsur- face signals within global SST syntheses would impact data-model comparison at the affected sites, and would also impact the amplitude of changes in global SST used in assessments of climate sensitivity.

Our analysis has highlighted the value of multi-proxy temperature reconstructions of the same samples from the same core, since this approach offers the best opportunity to explore the processes which could influence individual proxies, including changing mixed layer depth, seasonality, or upwelling intensity through time. A multi-proxy approach also allows the assessment of contributions by non-thermal factors to all of our proxy meth-
ods. Unfortunately, such an approach has rarely been applied at the sites we have examined here: the existing reso-
lution of single-site multi-proxy data prevented us from achieving our goal of evaluating in detail the processes
influencing differences and similarities in temperature records generated within and between sites. However, as
the majority of our records do record mean annual SSTs, with some caveats, we utilize all of the SST proxies we
have discussed here ($U_{37}^{K'}$, $TEX_{86}$, Mg/Ca) during our analysis and discussion in Section 4.

4. Assessing Trends and Variability in Ocean Temperatures, Ice Volume, and Sea Level Associated Across oNHG and iNHG

The transition from the late Pliocene to the early Pleistocene features ice-sheet expansion and a fall in pCO₂ (Figure 1; Section 1). An analysis of the LR04 δ¹⁸O_benthic stack (Lisiecki & Raymo, 2005) showed that there was a shift away from a relatively noisy relationship between δ¹⁸O_benthic and orbital forcing (a “stochastic” climate) ~2.8 Ma (Meyers & Hinnov, 2010). After 2.8 Ma, a closer relationship between δ¹⁸O_benthic and orbital forcing was identified (as a more predictable, or “deterministic” climate), and attributed to an increasing influence of larger ice sheets over climate feedbacks (Meyers & Hinnov, 2010). A change point in the same δ¹⁸O_benthic stack at 2.73 Ma (±0.1 Ma) was driven by an increase in the orbital obliquity signal (Ruggieri, 2013), and also highlights a shift in the climate forcing-response relationship. As the oceans are a critical part of heat and moisture transport through the climate system, interacting with the ice sheets in a variety of ways, we seek to identify their response to, or influence over, the iNHG.

The focus of this Section is to present our analysis of the high-resolution, well-dated records of SST, ice volume, and sea level in our PAGES-PlioVAR synthesis, which we refer to here as the “PlioVAR data sets” (Figure 2, Table 1). Our aim is to constrain the timings of events associated with the iNHG, and to understand the potential interactions between climate variability and longer-term climate evolution. As outlined in Section 1, we make a distinction between the oNHG (a broad window of time with evidence for slow ice-sheet growth) and iNHG (a more focused interval with a rapid increase IRD to the high latitude northern oceans). We first assess the statistical signatures of long-term trends (Section 4.1) and shorter-term variability with 41 kyr-obliquity frequency (Section 4.2) in the PlioVAR data sets. We then draw on complementary published data to evaluate and explore the processes driving the patterns of climate change we observe, by targeting the three research questions we posed in Section 1: Q2 What were the characteristics of the mPWP (~3.025–3.264 Ma) (Section 4.3); Q3 What were the characteristics of the iNHG? (Section 4.4); Q4 Did the amplitude or pacing of climate variability change with iNHG? (Section 4.5). Alongside our analysis of the PlioVAR data sets, we draw on additional marine sediment evidence to illustrate environmental changes which occurred onshore, particularly ice-sheet expansion (using IRD) (Andrews, 2000), and atmospheric processes via aeolian deposition of terrigenous materials which can reflect wind intensity, atmospheric circulation shifts, and sediment availability for transport for example, from aridity (Rea, 1993).

**4.1. PlioVAR Data-Sets: Long-Term Trends**

To test whether climate change during the iNHG had a globally synchronous onset, we employed change point
analysis on the SST, δ¹⁸O_benthic, δ¹⁸O_planktonic and sea-level data outlined in Section 2 (Figure 4, Table 3). The
change point analysis serves to test whether there have been statistically significant regime shifts, by identifying
the mean timing of the point which divides each record into two distinct intervals according to their mean state
or orbital-scale variability (for the methodology see Section 2.3). Since δ¹⁸O_benthic, δ¹⁸O_planktonic and sea-level data
should include a signature of global ice volume change, our null hypothesis was that a common change point
would detail the onset of globally significant ice-sheet growth consistent with the iNHG. In contrast, regional or
record-specific change points would reveal additional influences over δ¹⁸O including temperature and salinity,
which may be most strongly expressed in δ¹⁸O_planktonic data (Technical Box). Our SST records then give a contrast-
ing analysis of surface ocean properties, providing information about ocean circulation and temperature changes
which might be linked to other forcings (e.g. pCO₂ or gateway changes).

When considering the whole period of interest (3.3–2.4 Ma), a wide temporal range of change points is detected (Figure 4, Table 3). This is not a reflection of our narrow time window of interest (3.3–2.4 Ma), since analysis of SST and δ¹⁸O<sub>planktonic</sub> records from some of the same sites using a different change point approach and a broader time window also identified a wide spread of change points with no strong regional patterns (Kaboth-Bahr &
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Change point analysis of the four high time-resolution ocean proxies spanning the 3.3–2.4 Ma interval. (a) Summary map of change point results mapped by proxy. Symbol size represents the probability (large = high probability); (b) benthic δ¹⁸O stack (Ahn et al., 2017) with key marine isotope stages (MIS) indicated as per Figure 2; (c) Change points in sea-surface temperature (SST) records plotted by time and latitude; (d) Change point probabilities in the SST records plotted by time; (e) Change points in benthic and planktonic δ¹⁸O plotted by time and latitude; (f) Change point probabilities in benthic and planktonic δ¹⁸O records plotted by time.</figcaption></figure>
<table><caption>**Table 3**
Change Point Analysis for the 3.3–2.4 Ma Time Interval as Shown in Figure 4</caption><thead><tr><th rowspan="2">Site</th><th rowspan="2">Latitude (°N)</th><th rowspan="2">Ocean basin</th><th rowspan="2">Benthic δ<sup>18</sup>O</th><th colspan="3">Change point (Ma) (probability)</th></tr><tr><th>Planktonic δ<sup>18</sup>O</th><th>Planktonic Mg/Ca</th><th>U<sup>K'</sup><sub>37</sub> (except t = TEX<sub>86</sub>)</th></tr></thead><tbody><tr><td>ODP 907</td><td>69.24</td><td>Atlantic</td><td></td><td></td><td></td><td>3.17 (0.01)</td></tr><tr><td>ODP 982</td><td>57.52</td><td>Atlantic</td><td>2.62 (0.52)</td><td></td><td></td><td>2.99 (0.54)</td></tr><tr><td>U1417</td><td>56.96</td><td>Pacific</td><td></td><td></td><td></td><td>3.29 (0.02)</td></tr><tr><td>U1313</td><td>41.00</td><td>Atlantic</td><td>2.77 (0.28)</td><td>2.56 (0.57)</td><td>2.77 (0.74)</td><td>3.14 (0.32)</td></tr><tr><td>SPP</td><td>37.29</td><td>Medit.</td><td></td><td></td><td></td><td>2.89 (0.16)</td></tr><tr><td>ODP 978</td><td>36.23</td><td>Atlantic</td><td>2.88 (0.39)</td><td>3.27 (0.003)</td><td></td><td>3.29 (0.001)</td></tr><tr><td>ODP 1208</td><td>36.13</td><td>Pacific</td><td>3.00 (0.002)</td><td></td><td></td><td></td></tr><tr><td>ODP 1012</td><td>32.28</td><td>Pacific</td><td></td><td></td><td></td><td>3.07 (0.46)</td></tr><tr><td>ODP 1148</td><td>18.84</td><td>Pacific</td><td>2.48 (0.02)</td><td>Na</td><td></td><td></td></tr><tr><td>ODP 722</td><td>16.60</td><td>Indian</td><td>3.11 (0.86)</td><td></td><td></td><td>2.54 (0.43)</td></tr><tr><td>ODP 999</td><td>12.74</td><td>Atlantic</td><td>3.04 (0.23)</td><td>3.27 (0.005)</td><td>Na</td><td></td></tr><tr><td>ODP 1143</td><td>9.36</td><td>Pacific</td><td>2.52 (0.02)</td><td>Na</td><td>2.51 (0.77)</td><td>2.60 (0.77)</td></tr><tr><td>ODP 1241</td><td>5.84</td><td>Pacific</td><td>3.15 (0.37)</td><td></td><td>2.45 (0.005)</td><td></td></tr><tr><td>ODP 758</td><td>5.40</td><td>Indian</td><td>2.42 (0.37)</td><td></td><td></td><td></td></tr><tr><td>ODP 806</td><td>0.32</td><td>Pacific</td><td>2.80 (0.49)</td><td></td><td>2.41 (0.001)</td><td></td></tr><tr><td>ODP 849</td><td>0.18</td><td>Pacific</td><td>3.23 (0.001)</td><td></td><td></td><td></td></tr><tr><td>ODP 662</td><td>-1.40</td><td>Atlantic</td><td>2.73 (0.25)</td><td></td><td></td><td>3.09 (0.00)</td></tr><tr><td>ODP 846</td><td>-3.09</td><td>Pacific</td><td>3.11 (0.70)</td><td></td><td>Na</td><td>Na</td></tr><tr><td>ODP 214</td><td>-11.30</td><td>Indian</td><td>3.29 (0.001)</td><td>Na</td><td>2.42 (0.001)</td><td></td></tr><tr><td>U1463</td><td>-18.97</td><td>Indian</td><td></td><td></td><td></td><td>Na (t)</td></tr><tr><td>ODP 1082</td><td>-21.09</td><td>Atlantic</td><td>2.75 (0.06)</td><td></td><td></td><td>2.54 (0.77)</td></tr><tr><td>ODP 1264</td><td>-28.53</td><td>Atlantic</td><td>3.03 (0.1)</td><td></td><td></td><td></td></tr><tr><td>ODP 1267</td><td>-28.10</td><td>Atlantic</td><td>2.53 (0.04)</td><td></td><td></td><td></td></tr><tr><td>ODP 1087</td><td>-31.47</td><td>Atlantic</td><td></td><td></td><td></td><td>3.25 (0.14)</td></tr><tr><td>DSDP 593</td><td>-40.51</td><td>Pacific</td><td>2.73 (0.003)</td><td></td><td></td><td>Na</td></tr><tr><td>ODP 1088</td><td>-41.14</td><td>Atlantic</td><td>3.23 (0.003)</td><td></td><td></td><td></td></tr><tr><td>ODP 1090</td><td>-42.91</td><td>Atlantic</td><td></td><td></td><td></td><td>3.28 (0.01)</td></tr><tr><td>DSDP 594</td><td>-45.68</td><td>Pacific</td><td>2.77 (0.96)</td><td></td><td></td><td>3.0 (0.40)</td></tr><tr><td>ODP 1092</td><td>-46.41</td><td>Atlantic</td><td>3.15 (0.11)</td><td></td><td></td><td></td></tr><tr><td>DSDP 704</td><td>-46.88</td><td>Atlantic</td><td></td><td>2.77 (0.1)</td><td></td><td></td></tr></tbody></table>
<table><caption>**Table 3**
*Continued*</caption><thead><tr><th rowspan="2">Record</th><th>Globally integrated time series of sea-level using benthic δ<sup>18</sup>O</th></tr><tr><th>Change point (Ma)</th></tr></thead><tbody><tr><td>Berends et al. (2021)</td><td>Na</td></tr><tr><td>Miller et al. (2020)</td><td>Na</td></tr><tr><td>Rohling et al. (2014)</td><td>Na</td></tr></tbody></table>

*Note.* Site locations are shown in Figure 2. “Na” means data was analyzed but no statistical result was given. Blank spaces mean either data is not available or it is insufficient for analysis (cf. Table 1).

The probability that the change point is statistically significant is given in parentheses after each change point date.

Mudelsee, 2022). One $U_{37}^{K'}$ (DSDP Site 593), one TEX$_{86}$ (IODP Site U1463) and three $δ^{18}O_{planktonic}$ (ODP Sites 1143, 1148, and 214) records returned no change points. No change points were detected in the three sea-level time series (Table 3).

Our change point analysis challenges our long-held view that $δ^{18}O_{benthic}$ records represent a globally synchronous signal dominated by ice volume. Instead, local differences in bottom-water temperatures or seawater $δ^{18}O$ may be more important in the Pliocene and early Pleistocene records (Mudelsee & Raymo, 2005). We found a broad spread in $δ^{18}O_{benthic}$ change points, covering the full range from 3.29 to 2.42 Ma (Table 3, Figures 4a and 4e). This is consistent with a previous analysis of individual $δ^{18}O_{benthic}$ records using a statistical “ramp” to describe the start, end and amplitude of events linked to iNHG, which identified a broad window of changes to $δ^{18}O_{benthic}$ between 3.6 and 2.4 Ma (Mudelsee & Raymo, 2005). The earliest $δ^{18}O_{benthic}$ change points in the PlioVAR data sets occurred between ~3.3 and 3.0 Ma (our Interval P-1, Figure 1, Section 2.3), indicating that the climate system was already showing signs of change during the mPWP. The earliest mPWP $δ^{18}O_{benthic}$ change points developed from ~3.3 Ma, in the Southern Hemisphere, but these have low probabilities and suggest that a regime shift was unlikely to have occurred, or was weakly expressed in the data (Figures 4e and 4f). From ~3.2 Ma, $δ^{18}O_{benthic}$ change points with higher probabilities occurred in the northern low latitudes (<30°N) (Figures 4a, 4e, and 4f). Change points during the mPWP were also recorded in SSTs, first from ~3.3 Ma at sites influenced by the extent of subpolar surface water masses and upwelling systems and the Mediterranean Sea, but also with relatively low probability (Figures 4a, 4c, and 4d). From ~3.2 Ma, higher probability SST change points occurred in the mid- and high-latitudes of the Northern Hemisphere (Figures 4a, 4c, and 4d). After the mPWP, only two change points occurred between 3.0 and 2.8 Ma, both at ~2.9 Ma (Figures 4c and 4e) in the Mediterranean Sea ($δ^{18}O_{benthic}$, SST).

From 2.8 Ma onwards a wide latitudinal spread of $δ^{18}O_{benthic}$ change points occurred, first indicated by records influenced by mid-latitude climate change in both hemispheres. A cluster of change points between 2.8 and 2.7 Ma occurred in three sites which are bathed by the intermediate waters (~1,000–1,500 m water depth) formed in the surface ocean of the Subantarctic region (South-east Atlantic, ODP Site 1082, and South-west Pacific, DSDP Sites 593 and 594) (Figure 4e). These change points coincided with a $δ^{18}O_{planktonic}$ change point in the Subantarctic Atlantic (ODP Site 704), and occurred within a broader interval of southern mid- and high-latitude cooling as recorded by SSTs in the South-west Pacific and South-east Atlantic (Benguela) upwelling system at ~2.8–2.6 Ma (McClymont et al., 2016; Petrick et al., 2018). In the mid-latitude North Atlantic (IODP Site U1313) there are synchronous and high probability change points in *G. bulloides* Mg/Ca and $δ^{18}O_{benthic}$ at 2.77 Ma (noting that the $δ^{18}O_{planktonic}$ change point comes later, at 2.56 Ma).

Regional differences in the timing and character of climate evolution between low- and high-latitude signals across the iNHG have previously been identified (Kaboth-Bahr & Mudelsee, 2022; Ravelo et al., 2004). Our results confirm comparatively late change points in sites from the low-latitudes, given a final cluster of change points in $δ^{18}O_{benthic}$ from the South China Sea and Indian Ocean from 2.5 to 2.4 Ma (ODP Sites 1143, 1148, and ODP758), but acknowledging that for three of these sites the change point probabilities are very low (Figures 4a, 4d and 4e). A final cluster of low-latitude SST change points occurred between 2.6 and 2.4 Ma: all the change points in low-latitude Pacific and Indian Ocean Mg/Ca records occurred during this time interval (but with low probability), alongside high probability change points in $U_{37}^{K'}$-derived SSTs (Figure 4). The timing of these change points is broadly aligned with the pronounced glacial-interglacial cycles of MIS 96–100 (Figure 4).

## 4.2. PlioVAR Data-Sets: Orbital-Scale (Glacial-Interglacial) Variability

In the absence of synchronous change points (Section 4.1), we analyzed three broad time windows to compare orbital-scale SST and $δ^{18}O_{benthic}$ variability at the 41 kyr-obliquity frequency (Figure 1, Section 2.3): the sustained warmth of the mPWP (Interval P-1, 3.25–3.0 Ma), the transition toward intensified glaciations in the late Pliocene (Interval P-2, 3.0–2.7 Ma) and the main interval of intensified glaciations spanning the latest Pliocene to early Pleistocene (Interval P-3, 2.7–2.4 Ma) (Figure 5). Although the MIS M2 glacial stage is marked by cooling outside of the mPWP SST range at some sites (Table 4, Figure S6 in Supporting Information S1), its inclusion in the P-1 Interval does not impact the broad patterns in the data. To test whether iNHG was associated with climate cooling, we calculated the mean SSTs for each of the three time windows: in 19 ocean temperature records (17 sites), the mPWP (Interval P-1) was warmer than Interval P-3 (Table 5), but only 7 records (6 sites) showed cool-
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Violin plots (symmetrical vertical histograms) of glacial-interglacial sea-surface temperature (SST) amplitudes for each site, calculated as the maximum range within a moving 41-kyr window along the duration of each period. SSTs recorded by $U_{37}^{K'}$ ($u$ on the x-axis), Mg/Ca in planktonic foraminifera (m), and $TEX_{86}$ (t). Sites are ordered by latitude, from northernmost (left) to southernmost (right). Site locations are shown in Figure 2 and include Atlantic, Indian, and Pacific Ocean records. (a) a comparison of the evolution of mean SST variability across the three time intervals, extracted from panels (b–d); (b–d) frequency distributions of SST amplitudes within each time interval (colors). Box-whisker plots are superimposed on the frequency distribution data: means are joined between sites by the solid line. The 16th to 84th percentiles are indicated by the white vertical box, and the full range of data is shown by the shading for each site. For the frequency distributions which include the marine isotope stages M2 stage in panel (a) see Figure S5 in Supporting Information S1.</figcaption></figure>

ing which exceeded calibration error, and three sites recorded warming (Table 5). The largest cooling generally occurred in the mid-latitude sites, whereas small SST changes occurred in the low latitudes.

We explored the probability distribution of SST variability by latitude using “violin” plots for SST (Figure 5) and $\delta^{18}O_{benthic}$ (Figure 6). In each plot, the mean variability is delineated by the solid line joining the sites, whereas the amplitude of the variability is shown by the vertical extent of data per site, and the probability of each amplitude is shown by the width of the “violin.” For SST variability, a broad increase in the mean and the amplitude of SST variability
<caption>**Table 4**</caption>

Differences in Mean Sea-Surface Temperatures (SSTs) Between the MIS M2 Glacial Stage (3.282–3.308 Ma) and the KM5c Interglacial (3.197–3.214 Ma)

<table><thead><tr><th rowspan="2">Site</th><th rowspan="2">Lat. (°N)</th><th rowspan="2">Proxy</th><th colspan="3">KM5c</th><th colspan="3">M2</th><th rowspan="2">KM5c-M2 anomaly (°C)</th><th rowspan="2">KM5c-M2 uncertainty (°C)</th></tr><tr><th>Mean</th><th>2 s.d.</th><th>n</th><th>Mean</th><th>2 s.d.</th><th>n</th></tr></thead><tbody><tr><td>ODP 907</td><td>69.24</td><td>U<sup>K'</sup><sub>37</sub></td><td>Na</td><td>Na</td><td>0</td><td>9.4</td><td>2.1</td><td>22</td><td>Na</td><td>Na</td></tr><tr><td>ODP 982</td><td>57.52</td><td>U<sup>K'</sup><sub>37</sub></td><td>13.6</td><td>0.5</td><td>5</td><td>14.1</td><td>1.2</td><td>10</td><td>2.5</td><td>0.6</td></tr><tr><td>U1417</td><td>56.96</td><td>U<sup>K'</sup><sub>37</sub></td><td>9.5</td><td>Na</td><td>1</td><td>9.4</td><td>3.3</td><td>2</td><td>0.1</td><td>n/a</td></tr><tr><td>U1313</td><td>41.00</td><td>U<sup>K'</sup><sub>37</sub></td><td>20.4</td><td>0.9</td><td>3</td><td>16.7</td><td>0.4</td><td>7</td><td>3.7</td><td>0.7</td></tr><tr><td>U1313</td><td>41.00</td><td>Mg/Ca-b</td><td>22.6</td><td>1.0</td><td>2</td><td>21.1</td><td>0.9</td><td>8</td><td>1.5</td><td>1.0</td></tr><tr><td>SPP</td><td>37.29</td><td>U<sup>K'</sup><sub>37</sub></td><td>25.6</td><td>0.2</td><td>4</td><td>26.4</td><td>0.4</td><td>18</td><td>-0.8</td><td>0.2</td></tr><tr><td>ODP 978</td><td>36.23</td><td>U<sup>K'</sup><sub>37</sub></td><td>21.8</td><td>0.1</td><td>2</td><td>26.4</td><td>0.4</td><td>4</td><td>1.6</td><td>0.2</td></tr><tr><td>ODP 1208</td><td>36.13</td><td>U<sup>K'</sup><sub>37</sub></td><td>22.7</td><td>0.4</td><td>4</td><td>18.9</td><td>1.1</td><td>7</td><td>2.8</td><td>0.6</td></tr><tr><td>ODP 1012</td><td>32.28</td><td>U<sup>K'</sup><sub>37</sub></td><td>28.8</td><td>0.1</td><td>7</td><td>21.4</td><td>0.5</td><td>4</td><td>1.3</td><td>0.3</td></tr><tr><td>ODP 722</td><td>16.60</td><td>U<sup>K'</sup><sub>37</sub></td><td>30.5</td><td>0.2</td><td>9</td><td>27.7</td><td>0.5</td><td>9</td><td>1.2</td><td>0.2</td></tr><tr><td>ODP 999</td><td>12.74</td><td>Mg/Ca-s</td><td>24.0</td><td>0.3</td><td>3</td><td>23.0</td><td>0.9</td><td>19</td><td>1.0</td><td>0.4</td></tr><tr><td>ODP 1143</td><td>9.436</td><td>U<sup>K'</sup><sub>37</sub></td><td>28.0</td><td>0.0</td><td>3</td><td>30.4</td><td>0.3</td><td>3</td><td>0.2</td><td>0.2</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>Mg/Ca-r</td><td>29.9</td><td>0.3</td><td>6</td><td>29.6</td><td>Na</td><td>1</td><td>0.3</td><td>Na</td></tr><tr><td>ODP 1241</td><td>5.84</td><td>MgCa-s</td><td>27.7</td><td>0.1</td><td>2</td><td>27.2</td><td>0.4</td><td>5</td><td>0.5</td><td>0.2</td></tr><tr><td>ODP 806</td><td>0.32</td><td>Mg/Ca-s</td><td>28.6</td><td>0.2</td><td>3</td><td>28.7</td><td>1.0</td><td>4</td><td>-0.2</td><td>0.6</td></tr><tr><td>ODP 662</td><td>-1.39</td><td>U<sup>K'</sup><sub>37</sub></td><td>25.6</td><td>0.2</td><td>7</td><td>26.1</td><td>0.9</td><td>12</td><td>1.9</td><td>0.3</td></tr><tr><td>ODP 846</td><td>-3.09</td><td>U<sup>K'</sup><sub>37</sub></td><td>25.8</td><td>0.2</td><td>4</td><td>24.6</td><td>0.4</td><td>13</td><td>1.0</td><td>0.2</td></tr><tr><td>ODP 214</td><td>-11.30</td><td>Mg/Ca-r</td><td>26.5</td><td>Na</td><td>1</td><td>26.0</td><td>Na</td><td>1</td><td>0.5</td><td>Na</td></tr><tr><td>U1463</td><td>-18.97</td><td>TEX<sub>86</sub></td><td>28.5</td><td>0.6</td><td>4</td><td>28.1</td><td>0.8</td><td>4</td><td>0.3</td><td>0.7</td></tr><tr><td>ODP 1082</td><td>-21.09</td><td>U<sup>K'</sup><sub>37</sub></td><td>19.6</td><td>0.4</td><td>3</td><td>25.7</td><td>0.3</td><td>4</td><td>0.1</td><td>0.4</td></tr><tr><td>ODP 1087</td><td>-31.47</td><td>U<sup>K'</sup><sub>37</sub></td><td>19.6</td><td>0.5</td><td>7</td><td>17.7</td><td>1.0</td><td>5</td><td>1.9</td><td>0.6</td></tr><tr><td>DSDP 593</td><td>-40.51</td><td>U<sup>K'</sup><sub>37</sub></td><td>15.0</td><td>0.8</td><td>2</td><td>14.2</td><td>1.2</td><td>4</td><td>0.8</td><td>1.2</td></tr><tr><td>ODP 1090</td><td>-42.91</td><td>U<sup>K'</sup><sub>37</sub></td><td>13.2</td><td>0.8</td><td>4</td><td>9.1</td><td>0.3</td><td>4</td><td>4.1</td><td>0.5</td></tr><tr><td>DSDP 594</td><td>-45.68</td><td>U<sup>K'</sup><sub>37</sub></td><td>12.0</td><td>0.2</td><td>5</td><td>10.6</td><td>1.1</td><td>7</td><td>1.4</td><td>0.5</td></tr></tbody></table>

Note. SSTs were reconstructed using the $U^{K'}_{37}$ index, $TEX_{86}$ index, or Mg/Ca in three planktonic foraminifera species ($G. bulloides$, ($MgCa-b$), $T. sacculifer$ ($MgCa-s$), and $G. ruber$ ($MgCa-r$)). For each time interval, the mean, standard deviation ($2\sigma$) and the number of data points used ($n$) are shown. The KM5c-M2 anomaly is generated by subtracting the M2 mean from the KM5c mean, so that positive values equal a warming trend. The uncertainty in the KM5c-M2 SST difference is calculated as the sum of the uncertainties for each interval, calculated by dividing the $2\sigma$ range by the square root of the number of data points used in that time interval.

with increasing latitude was identified in all three time intervals (Figure 5), noting that the maximum latitude recorded for the Northern Hemisphere (ODP Site 907, 69°N) is greater than in the south (DSDP Site 594, 46°S). In contrast to this broad latitudinal pattern, the largest increase in mean variability through time occurred in the mid-latitude North Atlantic (IODP Site U1313, from Mg/Ca), where the range of variability was also large (as seen in both Mg/Ca in $G. bulloides$ and $U^{K'}_{37}$). Comparable latitudinal patterns of variability have been observed through the transition from the Last Glacial Maximum (~27–19 ka) into the Holocene interglacial (e.g. Rehfeld et al., 2018) and across multiple late Pleistocene glacial-interglacial cycles (Rohling et al., 2012). These patterns reflect higher sensitivity to global radiative forcing at high latitudes, alongside interactions between ice albedo (mid and high latitudes) and potentially hydrological feedbacks in the low latitudes, which can enhance mid-latitude SST variability (Rohling et al., 2012).

The Southern Hemisphere sites displayed small increases in the mean and range of orbital-scale SST variability after the mPWP. First, variability increased in the mid-latitude sites as early as the transition from Interval P-1 to P-2 (~3.0 Ma), especially in sites where the reconstructed SSTs are influenced by the position of the Subantarctic
<caption>**Table 5**</caption>

Quantifying the Mean Sea-Surface Temperature (SST) Differences Between the Three Time Intervals of Interest, Using SSTs Reconstructed by the $U_{37}^{K'}$ Index, $TEX_{86}$ Index, or Mg/Ca in Three Planktonic Foraminifera (G. bulloides, (MgCa-b), T. sacculifer (MgCa-s), and G. ruber (MgCa-r))

<table><thead><tr><th rowspan="2">Site</th><th rowspan="2">Lat.<br>(°N)</th><th rowspan="2">Proxy</th><th colspan="2">P-3 (2.7–2.4 Ma)</th><th colspan="2">P-2 (3.0–2.7 Ma)</th><th colspan="2">P-1 (3.3–3.0 Ma)</th><th rowspan="2">P3-P1<br>difference<br>(°C)</th><th rowspan="2">P3-P1<br>uncertainty<br>(°C)</th></tr><tr><th>Mean</th><th>2 s.d.</th><th>Mean</th><th>2 s.d.</th><th>Mean</th><th>2 s.d.</th></tr></thead><tbody><tr><td>ODP 907</td><td>69.24</td><td>$U_{37}^{K'}$</td><td>7.3</td><td>3.3</td><td>6.0</td><td>3.0</td><td>7.7</td><td>2.9</td><td>-0.4</td><td>0.9</td></tr><tr><td>ODP 982</td><td>57.52</td><td>$U_{37}^{K'}$</td><td>13.6</td><td>1.5</td><td>15.7</td><td>1.4</td><td>15.7</td><td>1.8</td><td>-2.1</td><td>0.3</td></tr><tr><td>U1417</td><td>56.96</td><td>$U_{37}^{K'}$</td><td>9.0</td><td>1.9</td><td>6.8</td><td>1.1</td><td>7.9</td><td>2.3</td><td>+1.1</td><td>0.9</td></tr><tr><td>U1313</td><td>41.00</td><td>$U_{37}^{K'}$</td><td>18.2</td><td>1.8</td><td>18.5</td><td>1.6</td><td>19.6</td><td>1.4</td><td>-1.4</td><td>0.3</td></tr><tr><td>U1313</td><td>41.00</td><td>Mg/Ca-b</td><td>18.6</td><td>2.7</td><td>19.5</td><td>3.0</td><td>21.2</td><td>1.3</td><td>-2.6</td><td>0.5</td></tr><tr><td>SPP</td><td>37.29</td><td>$U_{37}^{K'}$</td><td>25.3</td><td>1.0</td><td>26.2</td><td>1.1</td><td>26.1</td><td>0.8</td><td>-0.8</td><td>0.2</td></tr><tr><td>ODP 1012</td><td>32.28</td><td>$U_{37}^{K'}$</td><td>18.6</td><td>2.0</td><td>20.5</td><td>1.5</td><td>21.6</td><td>1.1</td><td>-3.0</td><td>0.3</td></tr><tr><td>ODP 722</td><td>16.60</td><td>$U_{37}^{K'}$</td><td>28.5</td><td>1.0</td><td>28.7</td><td>0.6</td><td>28.7</td><td>0.5</td><td>-0.2</td><td>0.1</td></tr><tr><td>ODP 999</td><td>12.74</td><td>MgCa-s</td><td>24.3</td><td>0.9</td><td>24.0</td><td>0.7</td><td>23.6</td><td>0.9</td><td>+0.7</td><td>0.3</td></tr><tr><td>ODP 1143</td><td>9.34</td><td>$U_{37}^{K'}$</td><td>29.6</td><td>0.5</td><td>30.0</td><td>0.4</td><td>30.3</td><td>0.3</td><td>-0.7</td><td>0.1</td></tr><tr><td>ODP 1143</td><td>9.36</td><td>Mg/Ca-r</td><td>29.7</td><td>0.5</td><td>29.6</td><td>0.6</td><td>29.8</td><td>0.5</td><td>-0.1</td><td>0.1</td></tr><tr><td>ODP 1241</td><td>5.84</td><td>Mg/Ca-s</td><td>27.0</td><td>0.6</td><td>27.5</td><td>0.5</td><td>27.6</td><td>0.4</td><td>-0.6</td><td>0.2</td></tr><tr><td>ODP 806</td><td>0.32</td><td>Mg/Ca-s</td><td>27.3</td><td>0.9</td><td>27.3</td><td>1.0</td><td>28.0</td><td>0.7</td><td>-0.7</td><td>0.3</td></tr><tr><td>ODP 662</td><td>-1.39</td><td>$U_{37}^{K'}$</td><td>27.3</td><td>1.0</td><td>27.6</td><td>0.8</td><td>27.8</td><td>0.9</td><td>-0.6</td><td>0.2</td></tr><tr><td>ODP 846</td><td>-3.09</td><td>$U_{37}^{K'}$</td><td>24.0</td><td>1.0</td><td>24.9</td><td>0.9</td><td>25.1</td><td>0.8</td><td>-1.0</td><td>0.2</td></tr><tr><td>ODP 214</td><td>-11.30</td><td>Mg/Ca-r</td><td>26.0</td><td>0.3</td><td>26.0</td><td>0.3</td><td>26.3</td><td>0.4</td><td>-0.3</td><td>0.2</td></tr><tr><td>U1463</td><td>-18.97</td><td>TEX<sub>86</sub></td><td>27.8</td><td>0.9</td><td>27.9</td><td>1.0</td><td>28.1</td><td>1.2</td><td>-0.4</td><td>0.4</td></tr><tr><td>ODP 1082</td><td>-21.09</td><td>$U_{37}^{K'}$</td><td>22.3</td><td>1.1</td><td>23.5</td><td>0.8</td><td>23.7</td><td>0.5</td><td>-1.4</td><td>0.2</td></tr><tr><td>ODP 1087</td><td>-31.47</td><td>$U_{37}^{K'}$</td><td>19.5</td><td>1.4</td><td>19.5</td><td>1.5</td><td>20.2</td><td>1.1</td><td>-0.6</td><td>0.4</td></tr><tr><td>DSDP 593</td><td>-40.51</td><td>$U_{37}^{K'}$</td><td>14.4</td><td>2.3</td><td>15.4</td><td>1.8</td><td>15.7</td><td>1.6</td><td>-1.4</td><td>0.8</td></tr><tr><td>ODP 1090</td><td>-42.91</td><td>$U_{37}^{K'}$</td><td>12.7</td><td>1.4</td><td>12.0</td><td>1.8</td><td>11.9</td><td>1.6</td><td>+0.8</td><td>0.4</td></tr><tr><td>DSDP 594</td><td>-45.68</td><td>$U_{37}^{K'}$</td><td>11.2</td><td>1.3</td><td>11.9</td><td>1.2</td><td>11.5</td><td>1.6</td><td>-0.2</td><td>0.6</td></tr></tbody></table>

Note. The P3-P1 anomaly is generated by subtracting the P-1 mean from the P3 mean, so that negative values equal early Pleistocene cooling relative to the Pliocene. The uncertainty in the P3-P1 SST difference is calculated as the sum of the uncertainties for each interval, calculated by dividing the $2\sigma$ range by the square root of the number of data points used in that time interval. Sites where the P3-P1 difference is less than the uncertainty are italicized. If the M2 is excluded from the P-3 interval (i.e., 3.25–3.0 Ma) the P-3 mean is changed by $\le 0.3^\circ C$, except at ODP Site 907 (0.9$^\circ C$ difference) and ODP Site 982 (0.4$^\circ C$ difference).

front in both the Pacific and Atlantic Oceans (e.g. ~41$^\circ$S, Southwest Pacific, DSDP Site 593) (e.g. McClymont et al., 2016). A second increase in SST variability occurred in the equatorial Pacific sites from Interval P-2 to P-3 (Figure 5; ODP Sites 806 and 846).

A more complex picture of SST variability emerged temporally and spatially in the Northern Hemisphere, which may be linked to a long-term expansion of the subpolar gyres (Martinez-Garcia et al., 2010; Sánchez-Montes et al., 2020) and intensification of subtropical circulation (Fedorov et al., 2015). Between Interval P-1 and P-2 (~3.0 Ma) there were increases in the mid- and high-latitude SST variability in several North Atlantic sites and in the Mediterranean Sea (Figure 5). Subsequently, moving from Interval P-2 to P-3 (~2.7 Ma) there was a large increase in SST variability in the mid-latitude Atlantic (IODP Site U1313), and in the North-east Pacific (Gulf of Alaska, IODP Site U1417; California Current, ODP Site 1012).

The same statistical analysis was applied to the $\delta^{18}O_{benthic}$ data (Figure 6). A systematic increase in both the amplitude and mean value of orbital-scale variability occurred from P-1 to P-3; Figure 6. This is clearly demonstrated in the LR04 $\delta^{18}O_{benthic}$ stack (Lisiecki & Raymo, 2005), which showed a small increase in mean variability at
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Violin plots (symmetrical vertical histograms) of glacial-interglacial δ¹⁸O_benthic amplitudes for each site, calculated as the maximum range within a moving 41-kyr window along the duration of each period, plotted using the same approach as for Figure 5. Colors refer to the water depth at the core site, except for the global LR04 stack (Lisiecki & Raymo, 2005) which is shown in black on the right.</figcaption></figure>

~3.0 Ma and a larger shift from ~2.7 Ma. However, the δ¹⁸O_benthic stack masks the diverse patterns of change in individual δ¹⁸O_benthic records, which are not easily grouped by latitude or water depth (Figure 6). For example, at 2,082 m water depth in the South Atlantic, ODP Site 1088 consistently showed the lowest mean variability in δ¹⁸O_benthic despite recording a wide range. In contrast, the deep ODP Site 849 (3,296 m water depth, equatorial Pacific) showed consistently low variability in both the mean and range throughout all of the time intervals we studied.

Although no change points were detected in the three sea-level records (Section 4.1), a consistent signal of increasing orbital-scale variability moving from Interval P-1 to Interval P-3 was recorded (Figure 7). The highest
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Violin plots (symmetrical vertical histograms) of glacial-interglacial sea-level amplitudes for each site, calculated as the maximum range within a moving 41-kyr window along the duration of each period, plotted using the same approach as for Figure 5. The records shown were generated using different approaches to extract sea-level change from the global $\delta^{18}O_{benthic}$ stack (Berends et al., 2021; Miller et al., 2020; Rohling et al., 2014) as outlined in Section 2.5.</figcaption></figure>

variability was consistently recorded in Miller2020 (~12–60 m), compared to ~10–50 m (Rohling2014) and 0–50 m (Berends2021). In the original sea-level time series, the Miller2020 reconstruction reached maximum variability from ~2.7 Ma, and all three sea-level reconstructions converged on a maximum variability at ~2.4 Ma (40–60 m) (Berends et al., 2021; Miller et al., 2020; Rohling et al., 2014).

### 4.3. What Were the Characteristics of the Mid-Piacenzian Warm Period (mPWP, ~3.264–3.025 Ma)? (Q2)

We have identified several SST and $\delta^{18}O_{benthic}$ change points which occurred within the mPWP, showing that this period included long-term changes in climate regime despite relative stability in the $pCO_2$ record (Figure 1). The earliest change points spanned both the MIS M2 glacial stage, which immediately preceded the mPWP, and the earliest interglacial of the mPWP (KM5c), that is, events were developing throughout our Interval P-1 (3.3–3.0 Ma). In this section, we compare our change point results with published records of ice-sheet size and configuration and complementary records of several key climate systems, to explore the nature of climate change occurring during both the M2 glacial stage and the mPWP.

The PlioVAR data sets reveal a wide range of SST responses to the MIS M2 glacial stage and the transition into interglacial KM5c (~3.25 Ma, Figures 1 and 5). The MIS M2 has previously been identified by a short-lived but pronounced increase in $\delta^{18}O_{benthic}$, argued to reflect the onset of more extensive glaciation in the Southern Hemisphere (e.g., McKay et al., 2012; Naish et al., 2009). There is also evidence for some ice-sheet growth in Greenland, Iceland and Svalbard during MIS M2 (De Schepper et al., 2014; Knies, Cabedo-Sanz, et al., 2014). However, despite MIS M2 being described as a possible “harbinger” of the onset of Pleistocene glacial cycles (Westerhold et al., 2020), there is a lack of widespread evidence for a large Northern Hemisphere ice-sheet expansion comparable to those of the Pleistocene (Kirby et al., 2020; Tan et al., 2017). Relative to the interglacial MIS KM5c, the PlioVAR data sets quantified the largest M2 cooling signals in SSTs from the mid- and high-latitudes of both hemispheres (2–4°C) (Table 4). Minimal cooling (<1°C) occurred in low-latitude SSTs (Table 4). This explains why the broad latitudinal increase in orbital-scale variability during the Interval P-1 was not altered if we included M2 in our analysis (Figure 5, Figure S5 in Supporting Information S1). Since only low probability change points in SST, $\delta^{18}O_{benthic}$ and $\delta^{18}O_{planktonic}$ occurred in the early part of Interval P-1, coinciding with MIS M2 and KM5c (3.3–3.2 Ma), we conclude that the changes in climate regime were weakly expressed. In contrast, from ~3.15 Ma change points with higher probabilities emerged, which were most significant in Northern Hemisphere SSTs and two low-latitude $\delta^{18}O_{benthic}$ records at ~3.1 Ma, during the second half of the mPWP.
Long-term shifts in regional circulation systems through the mPWP may be reflected in the wide latitudinal distribution of SST change points, which have then been captured locally with apparent asynchronicity in response to migrations of the positions of surface ocean fronts or upwelling systems through time. The SST change points which occur during the mPWP, and which precede $\delta^{18}O_{benthic}$ change points, suggest a rapid surface ocean response to forcing, either through direct radiative forcing or via locally sensitive shifts to surface ocean circulation systems. Early SST changes preceding $\delta^{18}O_{benthic}$ have also been observed in the mid- and high-latitude North Atlantic region over mid and late Pleistocene glacial-interglacial cycles in response to shifts in the positions of the Subpolar or Arctic Fronts (e.g. Alonso-Garcia et al., 2011; Barker et al., 2015; Hernández-Almeida et al., 2012, 2013; Mokeddem et al., 2014; Wright & Flower, 2002). In contrast, the later development of $\delta^{18}O_{benthic}$ change points will include the additional signal of the slower expansion of continental ice (e.g. Alonso-Garcia et al., 2011; Elderfield et al., 2012).

Although the western Mediterranean Sea SST and $\delta^{18}O_{benthic}$ change points had low probability (ODP Site 978, Figure 4), they aligned with reconstructed increases in west Mediterranean Sea surface salinity (~2 psu) and bottom water salinity (~1 psu; also ODP Site 978), which are proposed to have contributed to a strengthening of Mediterranean Outflow Water into the Atlantic Ocean by 3.3 Ma (Khélifi et al., 2014). The Mediterranean change points at ~3.3 Ma are also consistent with more restricted Indonesian Throughflow after ~3.5 Ma (Auer et al., 2019; De Vleeschouwer et al., 2022; Karas et al., 2009, 2011), which impacts Mediterranean SSTs and surface salinities by reducing African monsoon rainfall (Sarnthein et al., 2018). Increasing aridity in Eastern and Central Asia from ~3.6 Ma (Lu & Guo, 2014) or earlier (W. Zhang et al., 2018) have also been linked to a reduction in atmospheric water vapor as a result of a cooling climate.

Further evidence for longer-term cooling preceding, then continuing through, the mPWP comes from both high- and low-latitude regions. Expanding polar water mass extent in the Subantarctic South Atlantic from ~3.5 Ma (Martinez-Garcia et al., 2010) aligns with cooling and sea ice expansion within the mPWP in the Ross Sea (Riesselman & Dunbar, 2013). Enhanced dust deposition to the Subantarctic South Atlantic (ODP Site 1090) reflects intensification or equatorward displacement of southern hemisphere wind systems, likely causing cooling in the Benguela upwelling system (Southeast Atlantic: ODP Sites 1081, 1082, 1084, and 1087) (Martinez-Garcia et al., 2011; Petrick et al., 2018; Rosell-Melé et al., 2014). A more restricted Indonesian Throughflow would also have reduced heat transport through the Indian Ocean (Auer et al., 2019; De Vleeschouwer et al., 2022; Karas et al., 2009, 2011), and led to cooling in the southeast Atlantic Ocean (Karas et al., 2011; Rosell-Melé et al., 2014).

In the North Atlantic (ODP Site 982), SSTs cooled from ~3.6 Ma as Arctic sea-ice expansion was recorded at DSDP Site 610 (Karas et al., 2020 and references therein), and continued across the mPWP (ODP Site 982) (Lawrence et al., 2009). The slow growth in global ice volume starting as early as ~3.6 Ma (Mudelsee & Raymo, 2005) further highlights shifting high latitude climates before and through the mPWP. In several low-latitude sites, the M2 glacial stage marks the inflection point between a cooling trend starting from ~3.5 Ma and a subsequent warming through the mPWP: in the Atlantic (ODP Site 662), Pacific (ODP Site ODP846) (Herbert et al., 2010) and Indian (ODP Site 214) Oceans (Karas et al., 2009). The mPWP thus straddles an interval of long-term regional climate changes, which particularly connect the low latitudes and Southern Hemisphere sites.

Our analysis supports previous work indicating that the oNHG was a long-term transition, involving complex interactions between different climate system components and regions, potentially beginning before the M2 and extending through the mPWP (as early as ~3.6 Ma; Kleiven et al., 2002; Meyers & Hinnov, 2010; Mudelsee & Raymo, 2005). Orbital-scale $pCO_2$ reconstructions are lacking before M2, but available data provide values which were either comparable to, or slightly lower than, the mPWP (de la Vega et al., 2020; Pagani et al., 2010; Seki et al., 2010), indicating that $pCO_2$ forcing is unlikely to explain the complex climate trends observed here. Herbert et al. (2010) proposed that the long-term trends in low-latitude SSTs (which spanned ~300–500 kyr) could reflect a response to long-wavelength orbital forcing, favoring glacial-stage cooling during eccentricity minima. However, this mechanism does not explain signals of cooling from ~3.5 Ma, since they develop during an eccentricity maximum.

Alternatively, it has been proposed that pre-mPWP shifts in heat transport through the Southern Hemisphere could have resulted from a reduction in the Indonesian Throughflow (Auer et al., 2019; De Vleeschouwer et al., 2022). As the Indonesian Throughflow is influenced by regional sea level, both tectonic changes to the gateway config-
uration (Auer et al., 2019) and a fall in sea level in response to continental ice sheet growth (De Vleeschouwer et al., 2022) have been proposed to explain the reduction in the throughflow. In the high northern latitudes, the only likely major sources of IRD before ~2.7 Ma were Svalbard (Knies, Cabedo-Sanz, et al., 2014) and north-east Greenland (Bachem et al., 2017; Jansen et al., 2000). Hence, it is likely that a significant component of the increasing δ¹⁸O₄₉benthic across the wider iNHG interval (from ~3.6 Ma; Mudelsee & Raymo, 2005) includes Antarctic ice sheet expansion. IRD deposition in Southern Ocean cores also indicate that extensive Antarctic glaciation developed prior to the mPWP (Figure 1) (e.g. Cook et al., 2013; M. A. Hansen et al., 2015; Passchier, 2011; Patterson et al., 2014). Although warming and retreat of the West Antarctic ice sheet recorded in the Amundsen Sea starts from ~4.2 Ma and extends through the mPWP, more extensive West Antarctic glaciation from MIS M2 ~3.3 Ma is indicated by ice-proximal cores in the Ross Sea (McKay et al., 2012; Riesselman & Dunbar, 2013). The relatively early onset of Southern Hemisphere δ¹⁸O₄₉benthic change points through the mPWP may reflect ice sheet impacts on Southern Ocean circulation. However, the complexity of the Antarctic signals, showing both advance and retreat of different parts of the ice sheet through time, may have influenced the observed asynchronicity of the δ¹⁸O₄₉benthic change points. We return to these complexities in Section 4.5, where we assess the evidence for changing orbital-scale variability in our PlioVAR data sets and evidence for ice-volume change both during and after the mPWP.

## 4.4. What Were the Characteristics of the iNHG? (Q3)

In this Section, we draw on published records of changes in the continental ice sheets and their interactions with the oceans, to consider the broader context of the patterns of long-term changes we identified. We focus here on the long-term trends; in Section 4.5 we consider how orbital-scale climate variability may have changed with iNHG. First, we consider the evidence recorded in δ¹⁸O₄₉benthic and IRD, which suggest two intervals of ice-sheet growth: during the mPWP (i.e., before ~2.7 Ma) and an intensification of major ice-sheet growth between ~2.8 and 2.4 Ma. We then reflect on the evidence for changes occurring elsewhere within the wider climate system to consider whether mechanisms to explain these differences can be detected.

### 4.4.1. Growth of Continental Ice Prior to ~2.7 Ma (oNHG)

We have identified a wide range of δ¹⁸O₄₉benthic change points, which spanned our whole interval of interest (3.3–2.4 Ma) and with a wide range of probabilities. This result aligns with previous analysis of individual δ¹⁸O₄₉benthic data using a different approach (Mudelsee & Raymo, 2005). However, this broad interval of change contrasts with the change point analysis of the global LR04 δ¹⁸O₄₉benthic stack which placed the iNHG at 2.73 ± 0.1 Ma (Ruggieri, 2013), and with our analysis of sea-level variability which indicated an interval of globally significant ice growth between Intervals P-2 to P-3 (i.e., around ~2.7 Ma) (Figure 7). The temporal spread of δ¹⁸O₄₉benthic change points indicates that factors other than global ice volume have contributed to the site-specific δ¹⁸O₄₉benthic shifts. The properties of intermediate and deep waters are initially set by high-latitude surface ocean conditions, and we have identified several δ¹⁸O₄₉benthic during the mPWP which align with some early mid- and high-latitude SST change points (Figure 4). However, the early SST change points are diachronous and of varying probability. Therefore, although there may be an influence of high-latitude SST changes on some of the δ¹⁸O₄₉benthic signals (or vice versa), a direct connection is difficult to establish based on the data presented here.

Alternatively, the asynchronous advance and retreat of smaller ice sheets (Section 4.3) prior to ~2.7 Ma might account for the difficulties experienced in isolating signals of ice volume response to orbital forcing using δ¹⁸O₄₉benthic (Meyers & Hinnov, 2010). Evidence for growth of relatively small ice sheets includes reduced terrestrial organic matter and increased meltwater inputs to the Gulf of Alaska from ~3.1 Ma (IODP Site U1417), attributed to an expanding Cordilleran ice sheet which had not yet reached the sea to generate IRD (Sánchez-Montes et al., 2020, 2022). In the deep North-west Pacific (ODP Site 1208), by extracting the temperature influence from the δ¹⁸O₄₉benthic signal, a gradual increase in global ice volume was determined between ~3.15 and 2.73 Ma (Woodard et al., 2014). In the absence of extensive Northern Hemisphere Glaciation at this time, a growing Antarctic ice sheet was proposed to account for the >11 m sea-level fall by 2.73 Ma (Woodard et al., 2014). For the same time interval, long-term cooling of ~2°C in North Atlantic Deep Water (ODP Site 607) in turn suggests that there was cooling in the high latitude (northern) surface ocean where North Atlantic Deep Water is formed (Sosdian & Rosenthal, 2009). North Atlantic SST change points through 3.15–2.99 Ma provide support for decreasing mid- and high-latitude SSTs in this region (ODP Site 907, IODP Site U1313, ODP Site 982;
Figure 4). A challenge is to isolate whether these cooling patterns were conducive to ice-sheet growth, or whether they contributed a cooling signature to the individual records of $\delta^{18}O_{benthic}$, potentially masking evidence for the expansion of relatively small and dynamic ice sheets (Mudelsee & Raymo, 2005).

### 4.4.2. Intensification of Major Continental Ice Sheet Growth ~2.7–2.4 Ma (iNHG)

The $\delta^{18}O_{benthic}$ records reveal both an increased amplitude of orbital scale variability from ~2.7 Ma (P-3, Figure 6) and a clustering of change points ~2.8–2.7 Ma (Figure 4) which suggest a growing global impact of increasing ice volume in $\delta^{18}O_{benthic}$. The timing of these events broadly aligns with the iNHG (Section 1; Haug et al., 2005). Although the change points in $\delta^{18}O_{benthic}$ are clustered, they span several glacial-interglacial cycles rather than a relatively abrupt and synchronous event (Figure 4).

Asynchronicity in ice-sheet expansion from ~2.7 Ma is also evidenced by two separate increases in Northern Hemisphere IRD starting from ~2.7 Ma (Figure 1). It is overly simplistic to relate IRD abundance to glacial extent (Andrews, 2000), and we recognize that the early stages of ice-sheet expansion occur inland without any commensurate IRD deposition. However, abundant and widespread IRD deposition in high-latitude marine settings provides a powerful marker of continental-scale glaciation that has extended to sea level (Kleiven et al., 2002; Shackleton et al., 1984). The first dramatic increase in IRD abundance occurred during MIS G6 (~2.72 Ma; Figure 1) in sediments from all sectors of the Nordic Seas (Jansen et al., 2000; Jansen & Sjöholm, 1991; Knies et al., 2009), the subarctic northwest Pacific (Bailey et al., 2011) and as far south as ~53°N in the subpolar North Atlantic (Bailey et al., 2013; Blake-Mizen et al., 2019; Kleiven et al., 2002). The IRD increase starting from MIS G6 is supported by additional sedimentary evidence demonstrating that it reflected a synchronous expansion of glaciation in Greenland, Scandinavia, the British Isles, and parts of North America, and the advance of ice caps to their iceberg-calving margins (e.g. Chow et al., 1996; Kleiven et al., 2002; Lein et al., 2022; Thierens et al., 2012; Vanneste et al., 1995). A large increase in aeolian dust deposition in the mid-latitude North Atlantic during cold stages from ~2.72 Ma may also reflect strengthening of glaciogenic dust sources (Naafs et al., 2012) or a shift in regional wind fields and vegetation biomes (Lang et al., 2014) in North America following ice-sheet expansion.

The second large increase in IRD occurs from MIS G4 (~2.63 Ma, Figure 1) in sites recording ice expansion beyond the circum-Nordic Sea landmasses, suggesting a delayed expansion of ice sheets on North America (as originally proposed by Maslin et al., 1998). The onset of abundant IRD deposition in the Gulf of Alaska ~2.63 Ma (Sánchez-Montes et al., 2020) lies within the dating error of glacio-fluvial gravels marking the existence of the early and most extensive Cordilleran ice sheet (~2.64$^{+0.20}_{-0.18}$ Ma) (Hidy et al., 2013). In the North Atlantic, the first abundant IRD deposition in the center of the Last Glacial Maximum IRD belt (~50°N at IODP Site U1308) occurs during MIS G4 (Bailey et al., 2010) and later on the southernmost fringe (~41°N at IODP Site U1313) during MIS 100 at ~2.52 Ma (Bolton et al., 2010; Lang et al., 2014).

The diachronous expansion of IRD deposition in the subpolar North Atlantic between 2.72 and 2.52 Ma is not a product of iceberg survivability, because significant glacial excursions of Subarctic Front surface waters into the mid-latitude North Atlantic occurred earlier, during cold stages from MIS 104 (Bolton et al., 2018; Hennissen et al., 2014) and perhaps as early as MIS G6 (Naafs et al., 2010). Iceberg calving models for the Last Glacial Maximum suggest that the source of IRD to this region of the North Atlantic is the Gulf of St. Lawrence in midlatitude North America (Bigg & Wadley, 2001). The relatively late arrival of IRD deposition at IODP Site U1313 therefore likely reflects the first time during the iNHG that the proto-Laurentide Ice Sheet extended into the mid-latitudes of North America, consistent with other evidence for the timing of the onset of midlatitude glaciation here (Balco & Rovey, 2010; Shakun et al., 2016).

### 4.4.3. Climate Changes Accompanying the Intensification of Major Continental Ice Sheet Growth ~2.7–2.4 Ma (iNHG)

The range of change points and the diverse signatures of orbital-scale variability in $\delta^{18}O_{benthic}$ (Figure 6) are a reminder that factors other than ice volume are also required to explain the differences we observed between sites, which could include temperature or salinity-driven changes in seawater $\delta^{18}O$ in the regions of deep water formation. Only one SST record has a change point between 2.9 and 2.6 Ma (Figure 4; *G. bulloides* Mg/Ca from North Atlantic site IODP U1313). Yet, SST and dust records from both hemispheres detail signatures of temperature changes and intensification of westerly wind systems in the mid- and high-latitudes during this iNHG interval.
The reconstructed changes to North Pacific surface ocean conditions were conducive to preconditioning North-
ern Hemisphere ice-sheet growth, and impacted the configuration of deep ocean circulation. Warmer SSTs from
~2.7 Ma onwards (North Pacific sites ODP 882 and IODP U1417) would have provided a source of moisture to
enhance the growth of the North American ice sheets, alongside the cooler winter SSTs and associated sea-ice
expansion (Haug et al., 2005; Sánchez-Montes et al., 2020). These patterns reflect the development of the halo-
cline in the North Pacific, whereby a strong vertical gradient in salinity leads to isolation of deep waters from the
surface, encouraging elevated summer warming but also potentially enabling sequestration of atmospheric CO2
into the deep ocean (Haug et al., 2005). Model-data comparison suggests that as this strong halocline developed,
there was a cessation of the North Pacific Deep Water formation which occurred during the Pliocene (Burls
et al., 2017; Ford et al., 2022).

Through MIS G6 (~2.72 Ma), an increase in deep ocean connectivity between the Pacific and Atlantic Oceans is evidenced by the convergence of North Atlantic and North Pacific δ¹⁸O<sub>benthic</sub> values (Woodard et al., 2014). The pronounced cooling of North Atlantic Deep Water ~2.7 Ma (Site 607) marked the end of a long-term trend started in the mPWP (Sosdian & Rosenthal, 2009; Woodard et al., 2014). This deep water cooling is likely linked to the acceleration of surface ocean cooling in the North Atlantic at ~2.7 Ma (ODP Site 982, Lawrence et al., 2009) and an expansion of sea ice, which reached its maximum winter extent by ~2.6 Ma in the Fram Strait (Clotten et al., 2018; Knies, Cabedo-Sanz, et al., 2014).

While the Northern Hemisphere showed a general pattern of cooling over the ~2.7 Ma iNHG interval, the Subant-
arctic Atlantic was characterized by an interval of relative warmth (ODP Site 1090, Martinez-Garcia et al., 2010).
The intermediate waters which are formed in the Subantarctic surface ocean were also impacted. Change points
in δ¹⁸O<sub>benthic</sub> from sites bathed by Antarctic Intermediate Water at ~2.7 Ma (Figure 4) reflect both a shift in the
mean and a range of variability (Figure 6) in the Subantarctic surface ocean. In the South-west Pacific (DSDP Site
593), reconstructed Antarctic Intermediate Water temperatures confirm that a pronounced cooling occurred in the
zone of intermediate water formation during MIS G6, set within an extended interval of surface ocean cooling in
the mid-latitude South-west Pacific ~2.8–2.6 Ma (McClymont et al., 2016).

As the mid- and high latitudes cooled and climate became more variable, there was an intensification of wind-driven dust deposition in the mid-latitudes of both hemispheres and across multiple ocean basins from MIS G6 at ~2.73 Ma (Abell et al., 2021; Lang et al., 2014; Naafs et al., 2012). Both the marine sediment data and climate models suggest that globally synchronous equatorward shifts in the westerlies occurred in response to stronger atmospheric temperature gradients as the high latitudes cooled and ice volume increased (Abell et al., 2021; X. Li et al., 2015). In North America, MIS G6 marks a shift from a wetter-than-modern continental climate during the Pliocene to a more arid climate for the Pleistocene (Lang et al., 2014). A long-term increase in aridity in the lower latitudes has also been recorded from ~2.7 Ma and explained as a response to a globally cooler climate with reduced evaporation supplying monsoon rains (e.g. X. Li et al., 2015). In central Africa, a shift toward more arid conditions was followed by an increase in dust fluxes associated with the oNHG, suggested to reflect strengthening winds in response to the increasing latitudinal temperature gradients driven by ice-sheet growth (Crocker et al., 2022; de Menocal, 2004). Increasing aridity in North-east and South-east Africa also developed with iNHG (Liddy et al., 2016; A. K. Taylor et al., 2021), linked to long-term cooling in east Indian Ocean SSTs (A. K. Taylor et al., 2021).

Although we identified a wide range of climate changes aligned with the iNHG, we do not observe regime
shifts (change points) or changes to variability in the low-latitude regions ~2.7 Ma, consistent with a previous
proposal that the iNHG is a transition which largely reflects high-latitude climate change (Ravelo et al., 2004).
Even some mid-latitude sites show little long-term SST response to the iNHG ~2.7 Ma, for example, MIS G6
is relatively cool at North Atlantic Site IODP U1313, but otherwise forms part of a gradual long-term cooling
which began ~3.1 Ma (Naafs et al., 2012). In the Southern Hemisphere, long-term stability or gradual cool-
ing is observed in SSTs at several sites: in the Benguela upwelling system (Petrick et al., 2018; Rosell-Melé
et al., 2014), the subantarctic South Atlantic (Martinez-Garcia et al., 2010), and equatorial Pacific (Lawrence
et al., 2006). An impact of sea-level change (and orbital forcing) on the Indonesian Throughflow was demon-
strated by an intensification of the Leeuwin Current between 2.9 and 2.6 Ma, but this was strongly influenced by
ice sheet expansion (De Vleeschouwer et al., 2022) and so has a strong link to high-latitude changes.

The stronger signals of ice-sheet expansion and mid- and high-latitude climate changes with iNHG are consistent
with a fall in pCO₂ increasing the potential for ice sheets to expand into the mid-latitudes. A puzzle is to explain
why pCO₂ fell below mPWP values during MIS G10 (~2.8 Ma) whereas IRD increases are detected from MIS G6 (~2.72 Ma) (Figure 1). One possibility is that with the time resolution and spatial distribution of our existing records, we have not yet fully identified the expression of climate changes and ice-sheet evolution across these relatively short glacial intervals, including the pCO₂ reconstructions. Alternatively, the second lowering of pCO₂ during MIS G6 may have been required to facilitate the development of major ice-sheet growth (Figure 1), or a temperature threshold for positive ice-mass balance may finally have been crossed. Although much of our focus has been on the Northern Hemisphere ice sheets, several climate patterns associated with iNHG also suggest a climate system response to changes in the Antarctic ice sheet, including polar water expansion in both hemispheres (Haug et al., 2005; Martinez-Garcia et al., 2010) and shifting Southern Hemisphere circulation patterns (e.g. De Vleeschouwer et al., 2022), which might also have been connected to the monsoons, Mediterranean Outflow, and North Atlantic circulation (Sarnthein et al., 2018). In turn, it has been suggested that the increased Mediterranean Outflow from ~2.9 to 2.7 Ma might have exerted a negative impact on Northern Hemisphere ice sheet growth, delaying a response to North Atlantic cooling (Kaboth-Bahr et al., 2021).

Finally, our analysis identified a last clustering of change points focused largely on the low latitudes (<20°N/S) in both SST and δ¹⁸O_benthic data, which occurred after ~2.6 Ma and seem to be linked to the development of the first large glacial cycles (MIS 100-96; Figure 4e). The drivers of these final change points are unclear. Several are clustered around the Indo-Pacific warm pool (e.g. DSDP Site 214, ODP Site 806, ODP Site 1143, ODP Site 1148) or reflect changing zonal gradients in the Pacific Ocean (e.g. between ODP Site 806 and ODP Site 846). We acknowledge that caution is required since the Mg/Ca change points have very low probability (Figure 4), but these low-latitude signals may be hinting at a late response to iNHG or even an independent evolution of the Walker Circulation system. Previously, Kaboth-Bahr and Mudelsee (2022) identified high spatial and temporal variability in change points spanning the Pliocene-Pleistocene from sites influenced by Walker Circulation, and suggested that different ocean basins affected by this complex system evolved asynchronously from as early as 3.5 Ma. Our analysis supports this previous work, and shows that many of the high-latitude changes associated with iNHG occurred much earlier (before ~3.0 Ma), whereas the low latitudes tend to show a later long-term response aligned with the development of pronounced glacial maxima from MIS 100.

### 4.5. Did the Amplitude or Pacing of Climate Variability Change With iNHG? (Q4)

In this Section we investigate climate changes which occurred on orbital timescales (tens to hundreds of thousands of years), to consider whether iNHG impacted the relationships between climate forcings and feedbacks for different components of the climate system. We have identified that the mid- and high-latitudes had higher SST variability than the low-latitudes, which increased in association with the iNHG (Figure 5). Here, we explore the pacing of this variability, because it has been proposed that with the growth of larger Northern Hemisphere ice sheets, two impacts on the climate response to orbital forcing were experienced. First, that larger ice sheets were slower to respond to orbital forcing (Lawrence et al., 2009; Meyers & Hinnov, 2010), and second, that the ice sheets enhanced climate feedbacks related to changing ice albedo, meridional temperature gradients and winds (e.g. Herbert et al., 2015; Martinez-Boti et al., 2015). In particular we consider the relative importance of climate cycles related to obliquity pacing (~41 kyr) and precessional pacing (~19–23 kyr; see Section 1) by comparing our results to the published results of site- or proxy-specific time-series analysis.

There is considerable complexity in the pacing of late Pliocene climate, given that both obliquity and precession signals can be identified across a wide range of regions and climate proxies. Strong obliquity signals are evident in ocean temperatures (Dwyer et al., 1995; Lawrence et al., 2006, 2009; Naafs et al., 2012), records of both surface (e.g. Grundwell et al., 2008; Naafs et al., 2010) and deep (e.g. Lang et al., 2014) ocean circulation, and aeolian dust deposition (e.g. Crocker et al., 2022; Ding et al., 2002; Naafs et al., 2010). Ice-sheet responses to obliquity are indicated in records of δ¹⁸O_benthic (Lisiecki & Raymo, 2005), global sea-level (Naish & Wilson, 2009), and high northern latitude IRD (Shackleton et al., 1984), as well as sediment sequences detailing potential collapses of the West Antarctic ice sheet (Naish et al., 2009).

In contrast, precession signals dominate SST and δ¹⁸O_benthic records from the Mediterranean, reflecting a tight connection between local insolation and ventilation of bottom water masses (Herbert et al., 2015; Khélifi et al., 2014). Precession signals also dominate multiple low-latitude SST records (Herbert et al., 2010) and the wind- and river-driven delivery of African sediments to marine sediment sequences offshore, linked to changing African monsoon strength (Crocker et al., 2022; de Menocal, 2004). Significant precession signals are also found
in the high latitudes: in ice-marginal records from the marine-based East Antarctic ice sheet margin (Bertram et al., 2018; Cook et al., 2013; Patterson et al., 2014; Williams et al., 2010), and records of Antarctic Ice Sheet variability from the Scotia Sea (South Atlantic) (Reilly et al., 2021). Precession signals are also reflected strongly in the New Zealand sea level record, which is paced by Southern Hemisphere local insolation changes (Grant et al., 2019). The presence of both strong obliquity and precession signals indicates that even within the Antarctic Ice Sheet there were likely different regional catchment sensitivities to orbital forcing (Golledge et al., 2017), which might explain the challenges of isolating a strong link between orbital forcing and the global $\delta^{18}O_{benthic}$ stack prior to iNHG (Meyers & Hinnov, 2010).

From ~2.7 to 2.6 Ma (MIS G6-G4; Interval P-3), an increasing influence of obliquity pacing appears in several components of the climate system, and argued to reflect strengthening feedbacks associated with the expansion of continental ice sheets into the mid-latitudes (Section 4.3) (e.g. Herbert et al., 2015). In turn, these stronger connections may have strengthened the relationship between the LR04 $\delta^{18}O_{benthic}$ stack and orbital forcing, which then persisted until ~1.3 Ma (Meyers & Hinnov, 2010). For example, in the Mediterranean Sea, a shift away from precession-driven local insolation changes and the development of enhanced cooling during glacial stages with an obliquity period develops from ~2.8 Ma (Herbert et al., 2015). Obliquity continues to dominate dust deposition to the North Atlantic, but the relative timing shifts: having led $\delta^{18}O_{benthic}$ in the mPWP and late Pliocene, dust varies in-phase with $\delta^{18}O_{benthic}$ by ~2.8 Ma, suggesting a close link between dust supply and glaciation (Naafs et al., 2012).

However, we acknowledge that not all regions respond to the events marking the iNHG (MIS G6-G4) that is, the intensification of glacial cycles at 2.7–2.6 Ma (Section 4.4). Continued complexity in regional glacial-interglacial variability signals after 2.6 Ma is indicated by several records where precession, rather than obliquity, continues to dominate. These include Subantarctic Atlantic SSTs (ODP Site 1090; Martinez-Garcia et al., 2010), and meltwater discharges from a proto-Laurentide (North American) ice sheet into the Gulf of Mexico (Shakun et al., 2016). The precession-paced meltwater discharges are proposed to reflect a strong ice sheet mass balance response to local summer insolation (Shakun et al., 2016), yet the North American-sourced wind-driven dust deposition to the mid-latitude North Atlantic continued to be paced by obliquity and glaciation-driven sediment availability (Naafs et al., 2012). This apparent discrepancy in the observed pacing of ice-sheet dynamics during the iNHG might be explained by the fact that dust generation on North America during this time was dominated by non-glaciogenic processes (Lang et al., 2014). Alternatively, it may simply reflect that glaciogenic processes driving dust generation on North America during the iNHG operated at high latitudes under the influence of changes in axial tilt, whereas meltwater discharges from the mid-latitude margins were more strongly influenced by precession. This example demonstrates that complex and sometimes apparently contradictory evidence for climate pacing and forcing/response can be determined even in single regions, but further emphasizes the need for multi-proxy and multi-region assessments of climate evolution associated with the iNHG.

## 4.6. Summary and Future Outlook

### 4.6.1. SST and Ice Sheet Evolution Through the Onset and Intensification of Northern Hemisphere Glaciation (oNHG and iNHG)

In summary, we have found evidence for ice sheet growth, high latitude ocean changes, and low-high latitude teleconnections developing on both orbital and longer-term timescales through the mPWP and late Pliocene. The complex mixture of different orbital-scale variability being represented by different regions and climate variables, and their evolution through the mPWP and late Pliocene, may explain why the identification of common signals of climate and ice sheet evolution through this time window were not apparent in our change point analysis. The signals we observe also support the proposal that as well as this time interval being one with a stochastic character, the expansion of continental ice sheets with iNHG strengthened some low-high latitude climate feedbacks, and led to a more stable response to orbital forcing (Meyers & Hinnov, 2010). Yet even the Antarctic ice sheet showed variable regional responses to insolation, making the oNHG and iNHG challenging intervals to characterize, and perhaps an unexpectedly complex time window given the past research focus on the mPWP as an interval of equilibrium climate response to modern and near-future $pCO_2$.

We observed that some regions of the climate system, particularly the low latitudes, showed little response to the onset of major ice-sheet expansion ~2.7 Ma (iNHG). This pattern contrasts with the more pronounced
low-latitude cooling signatures which are observed during the Mid and Late Pleistocene glacial stages (e.g. Lawrence et al., 2006). As the impact of CO₂ forcing is highest in the low latitudes (Rohling et al., 2021), the relatively muted response to oNHG in these regions suggests that only a small decrease in pCO₂ occurred at this time. In contrast, the stronger expression of change in the mid- and high-latitudes is consistent with strengthening ice-albedo feedbacks which amplified the response to pCO₂ (Martinez-Boti et al., 2015; Rohling et al., 2012).

By combining these observations of various aspects of late Pliocene and early Pleistocene climate change, we have identified varying influences of ocean gateways, radiative forcing associated with pCO₂, ocean circulation and ice-sheet feedbacks. The early changes observed in sites which have connections to the Indonesian Through-flow, in the absence of any clear changes to pCO₂, suggest that changes to this ocean gateway may have had regional climate impacts prior to and during the mPWP. These events are accompanied by relatively localized and likely asynchronous patterns of growth in smaller ice masses, leaving a complex signature in both δ¹⁸Obenthic and SST signals. It is uncertain whether the ice sheets in both hemispheres were undergoing a long-term expansion during this interval, whereby shifts in ocean gateways and ocean circulation may have pre-conditioned a later expansion (iNHG), so that a small change in forcing led to a large ice-sheet response. Alternatively, the ice masses may have remained relatively small, in which case a much larger forcing would have been required to cause the larger expansion of continental ice defined by iNHG. Ice-proximal marine sediment records may offer some additional clues (see Section 4.6.2), but where ice margins remained on land there is a need to generate terrestrial data sets which can detect evidence for changes to air temperature, precipitation or ice margin position. Although some data exist (e.g. Hidy et al., 2013) it is generally sparsely distributed. An increase in terrestrial data density would also be valuable for data-model comparisons in warmer climates of the past (e.g. Haywood et al., 2020).

The later, and more abrupt, iNHG around 2.8–2.4 Ma has a stronger signal of northern hemisphere ice-sheet expansion as well as impacts on SSTs. The preceding fall in pCO₂ and latitudinal expressions of SST change are consistent with an ice-sheet response to reduced radiative forcing. The delayed response of the ice-sheets to the initial fall in pCO₂ at 2.8 Ma suggests either a role for slow ice-sheet feedbacks during growth, a need for further pCO₂ decline to have an impact, or some combination of the two. The asynchronous expansion of ice masses around, and then beyond, the Nordic Seas region suggests that there were different regional responses to these scenarios. Our analysis shows that once these larger ice masses grew, their impact on other parts of the climate system also grew, with some low latitude regions showing a response depending on the nature of the teleconnections at play.

### 4.6.2. Reflections and Future Work

In this work, we synthesized a large set of SST data and assessed their long-term trends and orbital-scale variability. The development of similar multi-proxy, spatially-distributed temperature data sets is also required to generate robust estimates of global mean surface temperature (GMST), which can be informative in data-model comparisons and calculations of climate sensitivity to changing radiative forcing (e.g. Inglis et al., 2020; Martinez-Boti et al., 2015). GMST estimates from MIS KM5c were derived by calculating a globally weighted mean temperature anomaly (see McClymont et al., 2020). However, this approach is sensitive to sampling density and does not account for unevenly distributed proxy data sets. Future studies can assess the influence of this approach by resampling the modern temperature field but with the sampling density of the Pliocene and Pleistocene. This has been assessed during other time intervals (e.g. the early Eocene) and indicates that this method is able to estimate GMST within 1.5°C (Inglis et al., 2020). Although δ¹⁸Obenthic values have also been used to estimate GMST at high-resolution during the Pliocene (e.g., J. Hansen et al., 2013), this approach is subject to several uncertainties (e.g., changes in ocean vertical stratification; discussed by Inglis et al., 2020). Furthermore, using δ¹⁸Obenthic values in this way assumes that they represent a global climate signal. Our study shows a wide range in δ¹⁸Obenthic change points starting from 3.3 Ma, suggesting that caution is required when using δ¹⁸Obenthic values to infer GMST, especially between glacial-interglacial cycles.

A limitation of our work was that to understand ice-sheet evolution, we relied on the expected global signature of ice volume from δ¹⁸Obenthic records, alongside published IRD evidence for ice sheets whose margins reach the sea. The diversity of the glacial-interglacial variability we identified between sites emphasizes the need to better understand both the global and local factors influencing the δ¹⁸Obenthic signal. More studies of the temperature-controlled Mg/Ca signal in benthic foraminifera would allow us to start to evaluate the varying environmental contributions to δ¹⁸Obenthic records. We identified some differences between the timing of changes recorded by δ¹⁸Obenthic and IRD (e.g. longer-term shifts in δ¹⁸Obenthic starting in the mid-Pliocene, a more
narrowly focussed shift in Northern Hemisphere IRD ~2.7–2.6 Ma). While several deep-sea ice-marginal records have provided continuous archives of ice-sheet variability, complexities exist in assessing the orbital pacing of ice sheets from these records. In part, this is because of the difficulty in separating out the specific processes controlling individual ice-sheet proxies (e.g. McKay et al., 2022). Detailed lithofacies analysis is required to ensure depositional processes are carefully considered, since otherwise there can be confusion even between studies undertaken at near-by sites (e.g. M. A. Hansen et al., 2015; Patterson et al., 2014).

For several regions of the ocean (e.g. the Indian Ocean, the Pacific sector of the Southern Ocean, the North Pacific Ocean), a paucity of published data also limits assessment of key systems for regulating regional and global climate, including monsoons and low-latitude hydroclimate, circulation systems connecting the low and high latitudes, and ocean-ice sheet interactions close to the ice sheets. Promising new high-resolution sequences recovered from the northern Indian Ocean by IODP will allow examination of the core monsoon region of the northern Bay of Bengal and Andaman Sea (Clemens et al., 2016). In the South-west Indian Ocean, IODP Expedition 361 has recovered sediment sequences offering new opportunities to reconstruct the transport of heat and salt from the low- to mid-latitudes and between the Indian and Atlantic Ocean, and to link them to onshore environmental changes in southern Africa (e.g. A. K. Taylor et al., 2021). A relative scarcity of data from the South Pacific sector of the Southern Ocean will soon be resolved as a result of IODP Expedition 383 (Lamy et al., 2021). Results from recent expeditions adjacent to the Antarctic Ice Sheet margin (IODP Expeditions 374, 379, 382) (Gohl et al., 2021; McKay et al., 2019; Weber et al., 2021), and scheduled expeditions adjacent to the Greenland margin (IODP Expedition 400) and close to Svalbard (Expedition 403) are also expected to shed new light on regional variability over high-latitude oceans and their proximal ice sheets. These expeditions will in future help to better inform the scientific community on the role of natural forcing and thresholds that may exist in warmer climates like the Pliocene Epoch compared to the relatively cool and more heavily glaciated Pleistocene Epoch. Reflecting on our analysis here, we also encourage the development of co-ordinated and strategic multi-proxy sampling strategies and generation of the highest quality age controls, to maximize the information which can be recovered and explored to better understand the iNHG and the warmth of the mPWP which immediately preceded it.

Finally, our focus on marine evidence to examine the iNHG neglects the rich archive of environmental information which can be returned from continental sequences, but whose age controls can be more challenging in the absence of time-continuous sedimentation or the global temporal framework of a $\delta^{18}O_{benthic}$ stack. Yet climate feedbacks linked to changes in vegetation and soils are important components of Earth system change (e.g. Feng et al., 2022). To maximize the potential for linking ocean, ice, and terrestrial information and gain a holistic understanding of the iNHG will require further analysis for continental environmental changes using both continental inputs to marine sediment sequences (cf., Crocker et al., 2022) and the generation and integration of continental sequences into a global framework (e.g. Feng et al., 2022). We also did not examine marine ecosystem responses to the environmental changes we identified, yet our interval of study includes evidence for diachronous extinction events (e.g. Brombacher et al., 2021; De Schepper & Head, 2008), and the mPWP offers an opportunity to consider how biota responded to a warmer climate with elevated $pCO_2$ (e.g. Todd et al., 2020).

**5. Conclusions**

By compiling, evaluating and exploring the high-resolution reconstructions of ocean temperatures, benthic and planktonic stable oxygen isotope signals ($\delta^{18}O_{benthic}$, $\delta^{18}O_{planktonic}$), and sea level we have identified complex signals of ice sheet and climate evolution between 3.3 and 2.4 Ma associated with the late Pliocene and earliest Pleistocene "onset and intensification of Northern Hemisphere Glaciation" (oNHG and iNHG). Here, we directly respond to the questions we posed in Section 1.

**Q1: Are there proxy-specific differences in the evolution of ocean temperatures?**

We did not identify any proxy-specific signatures where high-quality multi-proxy SST data were available. We confirmed that caution is required when considering whether data are recording mean annual SST: seasonality, depth of production, and how these two factors have evolved through time, as well as the possibility of non-thermal influences, are all important considerations when synthesizing and interpreting proxy data. The existing resolution of single-site multi-proxy data prevented us from fully addressing this question.
## Q2: What were the characteristics of the mPWP (~3.264–3.025 Ma)?

Immediately preceding the mPWP was the MIS M2 glacial stage, argued to reflect a short-lived increase in Southern Hemisphere glaciation and some minor ice-sheet expansion in the north. Our analysis shows that the M2-KM5c transition was marked by a reduction in global ice volume, and mid- and high-latitude ocean warming of 1–2 and 2–4°C, respectively. Low latitude ocean temperature response to the M2-KM5c transition was comparatively muted and within proxy calibration uncertainties (<1°C).

The mPWP was marked by warmer-than-Pre-Industrial climates, consistent with reduced global ice volume and elevated CO₂. It was also a time of climate evolution: sectors of the Antarctic ice sheet were expanding, and there is some evidence for growth in the smaller Northern Hemisphere ice sheets. There is a consistent pattern of increasing climate variability at the 41 kyr-obliquity scale moving from the low to the mid- and high-latitudes. Our analysis thus confirms that targeting a single interglacial for data-model comparison offers the best chance to minimize the influence of orbital-scale high-latitude ocean temperature variability for data-model comparison (Haywood et al., 2020; McClymont et al., 2020). In the absence of high-resolution *pCO₂* reconstructions before the mPWP, it is difficult to assess the role played by radiative forcing in the long-term trends. The observed patterns are consistent with ocean circulation responses to gateway-driven changes in Indonesian Throughflow and long-term intensification of meridional temperature gradients.

## Q3: What were the characteristics of the iNHG?

We have shown that the iNHG was not a globally synchronous event. We confirmed previous suggestions that long-term changes in climate and ice sheets began before and during the mPWP (a longer-term “oNHG”). Asynchronous transitions in δ¹⁸O_benthic records from individual sites confirm the importance of identifying and assessing the relative importance of local influences over δ¹⁸O_benthic (e.g. temperature, salinity) which are in turn linked to high-latitude surface ocean properties, alongside the more integrated signature of global ice volume.

Our review of existing data confirms that the iNHG was complex in terms of ice-sheet advances/retreats, in both timing and magnitude. A global expression of ice volume change is not indicative of how dynamic the individual ice masses were across the complete oNHG and iNHG. The pattern of high-latitude cooling proceeding first, followed by an expansion of ice towards the mid-latitudes in the northern hemisphere, is consistent with intensification of ice-albedo feedbacks which amplified the decline in CO₂ and enabled ice sheets to expand and survive further to the south. The lag between the fall in CO₂ during MIS G10 and the subsequent onset of IRD and glacial-stage intensification (Figure 1) may reflect the influence of changes in North Atlantic circulation (e.g. Kaboth-Bahr et al., 2021) or the time required to grow and sustain a larger continental ice mass, particularly in North America.

## Q4: Did the amplitude or pacing of climate variability change during iNHG?

Yes, but in a complex and non-uniform way. In almost all SST and δ¹⁸O_benthic records we examined, we observed an increase in the mean or range of variability at the 41 kyr-obliquity frequency from the mPWP through to the early Pleistocene, but there was strong regional scale variability. The low-latitude SST records show minimal or no shifts, whereas larger increases in SST variability emerge in the mid- and high-latitude Northern Hemisphere oceans. There was diversity in the signals of δ¹⁸O_benthic which may reflect evolving surface ocean conditions in the regions of intermediate and deep water formation, though no clear patterns could be isolated. In some proxy records, previous work has established switches from precession to obliquity in response to strengthened low-high latitude teleconnections as the Northern Hemisphere ice sheets expand, but not all records or regions express this pattern.

Our analysis has demonstrated the complex, non-uniform and globally asynchronous nature of climate changes associated with the oNHG and iNHG. The evidence available today indicates that shifting ocean gateways and ocean circulation changes may have pre-conditioned the later evolution of ice sheets with falling *pCO₂*, but there remains complexity in confirming the mechanism(s) or sequence of events which influenced the timing of ice sheet expansion and its expression in marine records. We recommend that future work targets multi-proxy data sets with high time-resolution, which can integrate our understanding of the environmental processes controlling and responding to the changes in SST and ice volume outlined here, which will continue to be critical to develop our understanding of the causes and impacts of past warm climates and the transitions between different climate regimes.
## Data Availability Statement

All original data-sets are available via the NOAA and PANGAEA online repositories and can be accessed through the citations and database DOIs given in Table 1, with full references provided below. The data used here, including any which were changed from the original publication due to differences in age control or proxy calibration, and the results of the analyses we performed, are available at https://doi.org/10.1594/PANGAEA.956158 (under moratorium). Data that we retrieved which was not on the NOAA or PANGAEA databases is also available at McClymont et al. (2023), https://doi.org/10.1594/PANGAEA.956158. The data we have used can also be accessed via our interactive web portal (https://pliovar.github.io/iNHG.html), including an option to download all data as a single file.

## Acknowledgments

We thank the Past Global Changes (PAGES) programme for their financial and logistical support of the working group on Pliocene Climate Variability over glacial-interglacial timescales (PlioVAR), and we thank all workshop participants and members of the PlioVAR steering committee for their discussions. Funding support has also been provided by the UK Natural Environment Research Council (PA, HLF [NE/N015045/1]), the Leverhulme Trust (Philip Leverhulme Prize, ELM), the Ministry of Science and Technology in Taiwan (SLH, Grants 110-2611-M-002-028-, 111-2611-M-002-027-), MICIN (MAG, the Grant PICTURE, PID2021-128322NB-100), a Royal Society Dorothy Hodgkin Fellowship (GNI, DHF/R1911178), the Research Council of Norway (BR, Grant 221712), the German Research Foundation (JG, project DePac-GR 3528/8-1), the Agencia Nacional de Investigación y Desarrollo, Anillo (CK, Grant ACT210046), the Universidad de Santiago de Chile (CK, Grant DICYT 092112KSSA), a UKRI Future Leadership Fellowship (BH, MR/S034293/1), the USA National Science Foundation (JA, Grant NSF-OCE-PRF #2126500) and the Agence Nationale de la Recherche in France (CTB, ANR-16-CE01-0004) and the IODP in the UK (PA) and France (CTB). This research used samples and/or data provided by the International Ocean Discovery Program (IODP), Ocean Drilling Program (ODP), and Deep Sea Drilling Project (DSDP). We also acknowledge PANGAEA and the NOAA Paleoclimatology databases for free access to published data. The work presented here resulted from the efforts of the late Antoni Rosell-Melé to generate a community effort to facilitate a better understanding of Pliocene and Pleistocene climate evolution: we thank him for all of his work and enthusiasm getting PlioVAR up and running, and for contributing to our discussions and assessments over the years that followed.

## References

Abell, J. T., Winckler, G., Anderson, R. F., & Herbert, T. D. (2021). Poleward sand weakened westerlies during Pliocene warmth. *Nature*, 589(7840), 70–75. https://doi.org/10.1038/s41586-020-03062-1

Ahn, S., Khider, D., Lisiecki, L. E., & Lawrence, C. E. (2017). A probabilistic Pliocene–Pleistocene stack of benthic $δ^{18}$O using a profile hidden Markov model. *Dynamics and Statistics of the Climate System*, 2(1), dxz002. https://doi.org/10.1093/climsys/dzx002

Alonso-Garcia, M., Sierro, F. J., & Flores, J. A. (2011). Arctic front shifts in the subpolar North Atlantic during the Mid-Pleistocene (800–400 ka) and their implications for ocean circulation. *Palaeogeography, Palaeoclimatology, Palaeoecology*, 311(3–4), 268–280. https://doi.org/10.1016/j.palaeo.2011.09.004

Anand, P., & Elderfield, H. (2005). Variability of Mg/Ca and Sr/Ca between and within the planktonic foraminifers Globigerina bulloides and Globorotalia truncatulinoidea. *Geochemistry, Geophysics, Geosystems*, 6(11), Q11D15. https://doi.org/10.1029/2004GC000811

Anand, P., Elderfield, H., & Conte, M. H. (2003). Calibration of Mg/Ca thermometry in planktonic foraminifera from a sediment trap time series. *Paleoceanography*, 18(2), 28–31. https://doi.org/10.1029/2002PA000846

Andersson, C., Warnke, D., Channell, J., Stoner, J., & Jansen, E. (2002). The mid-Pliocene (4.3–2.6 Ma) benthic stable isotope record of the Southern Ocean: ODP sites 1092 and 704, Meteor Rise, Palaeogeogr. *Palaeogeography, Palaeoclimatology, Palaeoecology*, 182(3–4), 165–181. https://doi.org/10.1016/S0031-0182(01)00494-1

Andrews, J. T. (2000). Icebergs and iceberg rated detritus (IRD) in the North Atlantic: Facts and assumptions. *Oceanography*, 13(3), 100–108. https://doi.org/10.5670/oceanog.2000.19

Auer, G., De Vleeschouwer, D., Smith, R. A., Bogus, K., Groeneveld, J., Grunert, P., et al. (2019). Timing and pacing of Indonesian through-flow restriction and its connection to late Pliocene climate shifts. *Paleoceanography and Paleoclimatology*, 34(4), 635–657. https://doi.org/10.1029/2018PA003512

Bachem, P. E., Risebrobakken, B., De Schepper, S., & McClymont, E. L. (2017). Highly variable Pliocene sea surface conditions in the Norwegian Sea. *Climate of the Past*, 13(9), 1153–1168. https://doi.org/10.5194/cp-13-1153-2017

Badger, M. P. S., Schmidt, D. N., Mackensen, A., & Pancost, R. D. (2013). High resolution alkenone palaeobarometry indicates stable $pCO_2$ during the Pliocene (3.3 to 2.8 Ma). *Proceedings of the Royal Society, A*, 371(2001), 20130094. https://doi.org/10.1098/rsta.2013.0094

Bailey, I., Bolton, C. T., DeConto, R. M., Pollard, D., Schiebel, R., & Wilson, P. A. (2010). A low threshold for North Atlantic ice rafting from "low-slung slippery" late Pliocene ice sheets. *Paleoceanography*, 25(1), PA1212. https://doi.org/10.1029/2009pa001736

Bailey, I., Hole, G. M., Foster, G. L., Wilson, P. A., Storey, C. D., Trueman, C. N., & Raymo, M. E. (2013). An alternative suggestion for the Pliocene onset of major northern hemisphere glaciation based on the geochemical provenance of North Atlantic Ocean ice-rafted debris. *Quaternary Science Reviews*, 75, 181–194. https://doi.org/10.1016/j.quascirev.2013.06.004

Bailey, I., Liu, Q., Swann, G. E., Jiang, Z., Sun, Y., Zhao, X., & Roberts, A. P. (2011). Iron fertilisation and biogeochemical cycles in the sub-Arctic northwest Pacific during the late Pliocene intensification of northern hemisphere glaciation. *Earth and Planetary Science Letters*, 307(3), 253–265. https://doi.org/10.1016/j.epsl.2011.05.029

Balco, G., & Rovey, C. W. (2010). Absolute chronology for major Pleistocene advances of the Laurentide Ice Sheet. *Geology*, 38(9), 795–798. https://doi.org/10.1130/g30946.1

Barker, S., Chen, J., Gong, X., Jonkers, L., Knorr, G., & Thornalley, D. (2015). Icebergs not the trigger for North Atlantic cold events. *Nature*, 520(7547), 333–336. https://doi.org/10.1038/nature14330

Barker, S., Greaves, M., & Elderfield, H. (2003). A study of cleaning procedures used for foraminiferal Mg/Ca paleothermometry. *Geochemistry, Geophysics, Geosystems*, 4(9), 8407. https://doi.org/10.1029/2003GC000559

Bé, A. W. H. (1980). Gametogenic calcification in a spinose planktonic foraminifer, Globigerinoides sacculifer (Brady). *Marine Micropaleontology*, 5, 283–310. https://doi.org/10.1016/0377-8398(80)90014-6

Bell, D. B., Jung, S. J. A., Kroon, D., Lourens, L. J., & Hodell, D. A. (2014). Local and regional trends in Plio-Pleistocene $δ^{18}$O records from benthic foraminifera. *Geochemistry, Geophysics, Geosystems*, 15(8), 3304–3321. https://doi.org/10.1002/2014GC005297

Berends, C. J., de Boer, B., & van de Wal, R. S. W. (2021). Reconstructing the evolution of ice sheets, sea level, and atmospheric $CO_2$ during the past 3.6 million years. *Climate of the Past*, 17(1), 361–377. https://doi.org/10.5194/cp-17-361-2021

Bertram, R. A., Wilson, D. J., van de Flierd, T., McKay, R. M., Patterson, M. O., Jimenez-Espejo, F. J., et al. (2018). Pliocene deglacial event timelines and the biogeochemical response offshore Wilkes Subglacial Basin, East Antarctica. *Earth and Planetary Science Letters*, 494, 109–116. https://doi.org/10.1016/j.epsl.2018.04.054

Bigg, G. R., & Wadley, M. R. (2001). The origin and flux of icebergs released into the Last Glacial Maximum Northern Hemisphere oceans: The impact of ice-sheet topography. *Journal of Quaternary Science*, 16(6), 565–573. https://doi.org/10.1002/jqs.628

Blake-Mizen, K., Hatfield, R. G., Stoner, J. S., Carlson, A. E., Xuan, C., Walczak, M., et al. (2019). Southern Greenland glaciation and Western boundary undercurrent evolution recorded on Eirik drift during the late Pliocene intensification of Northern Hemisphere glaciation. *Quaternary Science Reviews*, 209, 40–51. https://doi.org/10.1016/j.quascirev.2019.01.015

Bolton, C. T., Bailey, I., Friedrich, O., Tachikawa, K., Garidel-Thoron, T., Vidal, L., et al. (2018). North Atlantic midlatitude surface-circulation changes through the Plio-pleistocene intensification of Northern Hemisphere glaciation. *Paleoceanography and Paleoclimatology*, 33(11), 1186–1205. https://doi.org/10.1029/2018pa003412

Bolton, C. T., Wilson, P. A., Bailey, I., Friedrich, O., Beer, C. J., Becker, J., et al. (2010). Millennial-scale climate variability in the subpolar North Atlantic Ocean during the late Pliocene. *Paleoceanography*, 25(4), PA4218. https://doi.org/10.1029/2010pa001951
Bova, S., Rosenthal, Y., Liu, Z., Godad, S. P., & Yan, M. (2021). Seasonal origin of the thermal maxima at the Holocene and the last interglacial. *Nature*, **589**(7843), 548–553. https://doi.org/10.1038/s41586-020-03155-x

Brierley, C. M., Fedorov, A. V., Liu, Z., Herbert, T. D., Lawrence, K. T., & LaRiviere, J. P. (2009). Greatly expanded tropical warm pool and weakened Hadley circulation in the early Pliocene. *Science*, **323**(5922), 1714–1718. https://doi.org/10.1126/science.1167625

Brombacher, A., Wilson, P. A., Bailey, I., & Ezard, T. H. G. (2021). The dynamics of diachronous extinction associated with climatic deterioration near the Neogene/Quaternary boundary. *Paleoceanography and Paleoclimatology*, **36**(6), e2020PA004205. https://doi.org/10.1029/2020PA004205

Brown, S. J., & Elderfield, H. (1996). Variations in Mg/Ca and Sr/Ca ratios of planktonic foraminifera caused by postdepositional dissolution: Evidence of shallow Mg-dependent dissolution. *Paleoceanography*, **11**(5), 543–551. https://doi.org/10.1029/96pa01491

Burke, K. D., Williams, J. W., Chandler, M. A., Haywood, A. M., Lunt, D. J., & Otto-Bliesner, B. L. (2018). Pliocene and Eocene provide best analogs for near-future climates. *Proceedings of the National Academy of Sciences*, **115**(52), 13288–13293. https://doi.org/10.1073/pnas.1809600115

Burls, N. J., Fedorov, A. V., Sigman, D. M., Jaccard, S. L., Tiedemann, R., & Haug, G. H. (2017). Active Pacific meridional overturning circulation (PMOC) during the warm Pliocene. *Science Advances*, **3**(9), e1700156. https://doi.org/10.1126/sciadv.1700156

Caballero-Gill, R. P., Herbert, T. D., & Dowsett, H. J. (2019). 100-kyr paced climate change in the Pliocene warm period, Southwest Pacific. *Paleoceanography and Paleoclimatology*, **34**(4), 524–545. https://doi.org/10.1029/2018pa003496

Capron, E., Govin, A., Stone, E. J., Masson-Delmotte, V., Mulitza, S., Otto-Bliesner, B., et al. (2014). Temporal and spatial structure of multimillennial temperature changes at high latitudes during the Last Interglacial. *Quaternary Science Reviews*, **103**, 116–133. https://doi.org/10.1016/j.quascirev.2014.08.018

Chen, J., Farrell, J. W., Murray, D. W., & Prell, W. L. (1995). Timescale and paleoceanographic implications of a 3.6 m.y. oxygen isotope record from the northeast Indian Ocean (Ocean Drilling Program Site 758). *Paleoceanography*, **10**(1), 21–48. https://doi.org/10.1029/94PA02290

Cheng, X., Zhao, Q., Wang, J., Jian, Z., Xia, P., Huang, B., et al. (2004). Data report: Stable isotopes from sites 1147 and 1148. In W. L. Prell, P. Blum, E. M. Arnold, C. Bühring, M. P. Chen, & S. Clemens, et al. (Eds.), *Proceedings of the ocean drilling program, scientific results* (Vol. 184, pp. 1–12). Ocean Drilling Program. https://doi.org/10.2973/odp.proc.str.184.223.2004

Chow, N., Morad, S., & Al-Aasm, I. S. (Eds.). (1996). *Origin of authigenic carbonates in Eocene to Quaternary sediments from the Arctic ocean and Norwegian Greenland sea* (pp. 415–434). Ocean Drilling Program.

Clemens, S. C., Kuhnt, W., & LeVay, L. J. (2016). Indian monsoon rainfall. In *Proceedings of the International Ocean Discovery program* (Vol. 353). International Ocean Discovery Program. https://doi.org/10.14379/iodp.proc.353.2016

Clemens, S. C., Murray, D. W., & Prell, W. L. (1996). Nonstationary phase of the Plio-pleistocene Asian monsoon. *Science*, **274**(5289), 943–948. https://doi.org/10.1126/science.274.5289.943

Clotten, C., Stein, R., Fahl, K., & De Schepper, S. (2018). Seasonal sea ice cover during the warm Pliocene: Evidence from the Iceland sea (ODP site 907). *Earth and Planetary Science Letters*, **481**, 61–72. https://doi.org/10.1016/j.epsl.2017.10.011

Conte, M. H., Sicre, M. A., Rühlemann, C., Weber, J. C., Schulte, S., Schulz-Bull, D., & Blanz, T. (2006). Global temperature calibration of the alkenone unsaturation index (UK') in surface waters and comparison with surface sediments. *Geochemistry, Geophysics, Geosystems*, **7**(2), Q02005. https://doi.org/10.01029/02005GC001054

Cook, C. P., van de Flierdt, T., Williams, T., Hemming, S. R., Iwai, M., Kobayashi, M., et al. (2013). Dynamic behaviour of the East Antarctic ice sheet during Pliocene warmth. *Nature Geoscience*, **6**(9), 765–769. https://doi.org/10.1038/ngeo1889

Craig, H., & Gordon, L. I. (1965). Deuterium and oxygen 18 variations in the ocean and the marine atmosphere. In E. Tongiogi & V. Lishi (Eds.), *Proceedings of a conference on stable isotopes in oceanographic studies and Paleotemperatures, Spoleto, Italy* (pp. 9–130).

Crocker, A. J., Naafs, B. D. A., Westerhold, T., James, R. H., Cooper, M. J., Röhl, U., et al. (2022). Astronomically controlled aridity in the Sahara since at least 11 million years ago. *Nature Geoscience*, **15**(8), 671–676. https://doi.org/10.1038/s41561-022-00990-7

Crundwell, M., Scott, G., Naish, T., & Carter, L. (2008). Glacial-interglacial ocean climate variability from planktonic foraminifera during the Mid-Pleistocene transition in the temperate Southwest Pacific, ODP Site 1123. *Palaeogeography, Palaeoclimatology, Palaeoecology*, **260**(1–2), 202–229. https://doi.org/10.1016/j.palaeo.2007.08.023

de la Vega, E., Chalk, T. B., Wilson, P. A., Bysani, R. P., & Foster, G. L. (2020). Atmospheric CO₂ during the mid-Piacenzian warm period and the M2 glaciation. *Scientific Reports*, **10**(1), 11002. https://doi.org/10.1038/s41598-020-67154-8

De Conto, R. M., Pollard, D., Wilson, P. A., Pälike, H., Lear, C. H., & Pagani, M. (2008). Thresholds for Cenozoic bipolar glaciation. *Nature*, **455**(7213), 652–656. https://doi.org/10.1038/nature07337

Dekens, P. S., Lea, D. W., Pak, D. K., & Spero, H. J. (2002). Core top calibration of Mg/Ca in tropical foraminifera: Refining paleotemperature estimation. *Geochemistry, Geophysics, Geosystems*, **3**(4), 1022–1029. https://doi.org/10.1029/2001GC000200

Dekens, P. S., Ravelo, A. C., & McCarthy, M. D. (2007). Warm upwelling regions in the Pliocene warm period. *Paleoceanography*, **22**(3), PA3211. https://doi.org/10.1029/2006PA001394

de Menocal, P. B. (2004). African climate change and faunal evolution during the Pliocene-Pleistocene. *Earth and Planetary Science Letters*, **220**(1–2), 3–24. https://doi.org/10.1016/s0012-821x(04)00003-2

De Schepper, S., Gibbard, P. L., Salzmann, U., & Ehlers, J. (2014). A global synthesis of the marine and terrestrial evidence for glaciation during the Pliocene Epoch. *Earth-Science Reviews*, **135**, 83–102. https://doi.org/10.1016/j.earscirev.2014.04.003

De Schepper, S., Groeneveld, J., Naafs, B. D. A., Van Renterghem, C., Hennissen, J., Head, M. J., et al. (2013). Northern Hemisphere glaciation during the globally warm early late Pliocene. *PLoS One*, **8**(12), e81508. https://doi.org/10.81371/journal.pone.0081508

De Schepper, S., & Head, M. J. (2008). Age calibration of dinoflagellate cyst and acritarch events in the Pliocene-Pleistocene of the eastern North Atlantic (DSDP Hole 610A). *Stratigraphy*, **5**(2), 137–161.

De Vleeschouwer, D., Peral, M., Marchegiano, M., Füllberg, A., Meinicke, N., Pälike, H., et al. (2022). Plio-Pleistocene Perth Basin water temperatures and Leeuwin Current dynamics (Indian Ocean) derived from oxygen and clumped-isotope paleothermometry. *Climate of the Past*, **18**(5), 1231–1253. https://doi.org/10.5194/cp-18-1231-2022

Ding, Z. L., Derbyshire, E., Yang, S. L., Yu, Z. W., Xiong, S. F., & Liu, T. S. (2002). Stacked 2.6-Ma grain size record from the Chinese loess based on five sections and correlation with the deep-sea δ¹⁸O record. *Paleoceanography*, **17**(3), 1033–5-21. https://doi.org/10.1029/2001PA000725

Dowsett, H., Dolan, A., Rowley, D., Moucha, R., Forte, A. M., Mitrovica, J. X., et al. (2016). The PRISM4 (mid-Piacenzian) paleoenvironmental reconstruction. *Climate of the Past*, **12**(7), 1519–1538. https://doi.org/10.5194/cp-12-1519-2016

Dowsett, H. J., Robinson, M. M., Haywood, A. M., Hill, D. J., Dolan, A. M., Stoll, D. K., et al. (2012). Assessing confidence in Pliocene sea surface temperatures to evaluate predictive models. *Nature Climate Change*, **2**(5), 365–371. https://doi.org/10.1038/nclimate1455

Dumitru, O.-A., Austermann, J., Polyak, V. J., Fornós, J. J., Asmerom, Y., Ginés, J., et al. (2019). Constraints on global mean sea level during Pliocene warmth. *Nature*, **574**(7777), 233–236. https://doi.org/10.1038/s41586-019-1543-2
Dupont, L. M., Donner, B., Vidal, L., Pérez, E. M., & Wefer, G. (2005). Linking desert evolution and coastal upwelling: Pliocene climate change in Namibia. *Geology*, *33*(6), 461–464. https://doi.org/10.1130/G21401.1

Dwyer, G. S., Cronin, T. M., Baker, P. A., Raymo, M. E., Buzas, J. S., & Corrège, T. (1995). North Atlantic deepwater temperature change during late pleistocene and late Quaternary climatic cycles. *Science*, *270*, 1247–1351.

Elderfield, H., Ferretti, P., Greaves, M., Crowhurst, S., McCave, I. N., Hodell, D., & Piotrowski, A. M. (2012). Evolution of ocean temperature and ice volume through the mid-Pleistocene climate transition. *Science*, *337*(6095), 704–709. https://doi.org/10.1126/science.1221294

Elderfield, H., & Ganssen, G. (2000). Past temperature and δ¹⁸O of surface ocean waters inferred from foraminiferal Mg/Ca ratios. *Nature*, *405*(6785), 442–445. https://doi.org/10.1038/35013033

Epstein, S., Buchsbaum, R., Lowenstam, H. A., & Urey, H. C. (1953). Revised carbonate-water isotopic temperature scale. *Geological Society of America Bulletin*, *64*(11), 1315–1326. https://doi.org/10.1130/0016-7606(1953)64[1315:rcits]2.0.co;2

Etourneau, J., Martinez, P., Blanz, T., & Schneider, R. (2009). Pliocene-Pleistocene variability of upwelling activity, productivity, and nutrient cycling in the Benguela region. *Geology*, *37*(10), 871–874. https://doi.org/10.1130/g25733a.1

Evans, D., Brierley, C., Raymo, M. E., Erez, J., & Müller, W. (2016). Planktic foraminifera shell chemistry response to seawater chemistry: Pliocene–Pleistocene seawater Mg/Ca, temperature and sea level change. *Earth and Planetary Science Letters*, *438*, 139–148. https://doi.org/10.1016/j.epsl.2016.01.013

Fairbanks, R. G., Sverdlove, M., Free, R., Wiebe, P. H., & Bé, A. W. H. (1982). Vertical distribution and isotopic fractionation of living planktonic foraminifera from the Panama Basin. *Nature*, *298*(5877), 841–844. https://doi.org/10.1038/298841a0

Fedorov, A. V., Burls, N. J., Lawrence, K. T., & Peterson, L. C. (2015). Tightly linked zonal and meridional sea surface temperature gradients over the past five million years. *Nature Geoscience*, *8*(12), 975–980. https://doi.org/10.1038/ngeo2577

Feng, R., Bhattacharya, T., Otto-Bliesner, B. L., Brady, E. C., Haywood, A. M., Tindall, J. C., et al. (2022). Past terrestrial hydroclimate sensitivity controlled by Earth system feedbacks. *Nature Communications*, *13*(1), 1306. https://doi.org/10.1038/s41467-022-28814-7

Fischer, H., Meissner, K. J., Mix, A. C., Abram, N. J., Austermann, J., Brovkin, V., et al. (2018). Palaeoclimate constraints on the impact of 2°C anthropogenic warming and beyond. *Nature Geoscience*, *11*(7), 474–485. https://doi.org/10.1038/s41561-018-0146-0

Ford, H. L., Burls, N. J., Jacobs, P., Jahn, A., Caballero-Gill, R. P., Hodell, D. A., & Fedorov, A. V. (2022). Sustained mid-Pliocene warmth led to deep water formation in the North Pacific. *Nature Geoscience*, *15*(8), 658–663. https://doi.org/10.1038/s41561-022-00978-3

Ford, H. L., Ravelo, A. C., Dekens, P. S., LaRiviere, J. P., & Wara, M. W. (2015). The evolution of the equatorial thermocline and the early Pliocene El Padre mean state. *Geophysical Research Letters*, *42*(12), 4878–4887. https://doi.org/10.1002/2015GL064215

Friedrich, O., Wilson, P. A., Bolton, C. T., Beer, C. J., & Schiebel, R. (2013). Late Pliocene to early Pleistocene changes in the North Atlantic Current and suborbital-scale sea-surface temperature variability. *Paleoceanography*, *28*(2), 274–282. https://doi.org/10.1002/palo.20029

García-Gallardo, Á., Grunert, P., & Piller, W. E. (2018). Variations in Mediterranean–Atlantic exchange across the late Pliocene climate transition. *Climate of the Past*, *14*(3), 339–350. https://doi.org/10.5194/cp-14-339-2018

Gasson, E., Lunt, D. J., DeConto, R., Goldner, A., Heinemann, M., Huber, M., et al. (2014). Uncertainties in the modelled CO₂ threshold for Antarctic glaciation. *Climate of the Past*, *10*(2), 451–466. https://doi.org/10.5194/cp-10-451-2014

Gohl, K., Uenzelmann-Neben, G., Gille-Petzold, J., Hillenbrand, C., Klages, J. P., Bohaty, S. M., et al. (2021). Evidence for a highly dynamic West Antarctic ice sheet during the Pliocene. *Geophysical Research Letters*, *48*(14), e2021GL093103. https://doi.org/10.1029/2021GL093103

Golledge, N. R., Levy, R. H., McKay, R. M., & Naish, T. R. (2017). East Antarctic ice sheet most vulnerable to Weddell Sea warming. *Geophysical Research Letters*, *44*(5), 2343–2351. https://doi.org/10.1002/2016GL072422

Gradstein, F. M., Ogg, J. G., Schmitz, M., & Ogg, G. (Eds.). (2012). *The geological time scale 2012*. Elsevier.

Grant, G. R., & Naish, T. R. (2021). Pliocene sea level revisited: Is there more than meets the eye? *Past Global Changes Magazine*, *29*(1), 34–35. https://doi.org/10.22498/pages.29.1.34

Grant, G. R., Naish, T. R., Dunbar, G. B., Stocchi, P., Kominz, M. A., Kamp, P. J. J., et al. (2019). The amplitude and origin of sea-level variability during the Pliocene epoch. *Nature*, *574*(7777), 237–241. https://doi.org/10.1038/s41586-019-1619-2

Gray, W. R., & Evans, D. (2019). Nonthermal influences on Mg/Ca in planktonic foraminifera: A review of culture studies and application to the last glacial maximum. *Paleoceanography and Paleoclimatology*, *34*(3), 306–315. https://doi.org/10.1029/2018pa003517

Groeneveld, J., Hathorne, E., Steinke, S., DeBey, H., Mackensen, A., & Tiedemann, R. (2014). Glacial induced closure of the Panamanian gateway during marine isotope stages (MIS) 96-101 (2.5 Ma). *Earth and Planetary Science Letters*, *404*, 296–306. https://doi.org/10.1016/j.epsl.2014.08.007

Groeneveld, J., & Tiedemann, R. (2005). Mg/Ca and sea surface temperature data for Site 202-1241. *PANGAEA*. https://doi.org/10.1594/PANGAEA.315654

Gulev, S. K., Thorne, P. W., Ahn, J., Dentener, F. J., Domingues, C. M., Gerland, S., et al. (2021). Changing state of the climate system. In V. Masson-Delmotte, P. Zhai, A. Pirani, S. L. Connors, C. Péan, S. Berger, et al. (Eds.), *Climate change 2021: The physical science basis. Contribution of working group I to the sixth assessment report of the intergovernmental panel on climate change* (pp. 287–422). Cambridge University Press.

Hansen, J., Sato, M., Russell, G., & Kharecha, P. (2013). Climate sensitivity, sea level and atmospheric carbon dioxide. *Philosophical Transactions of the Royal Society A: Mathematical, Physical & Engineering Sciences*, *371*(2001), 20120294. https://doi.org/10.1098/rsta.2012.0294

Hansen, M. A., Passchier, S., Khim, B., Song, B., & Williams, T. (2015). Threshold behavior of a marine-based sector of the East Antarctic Ice Sheet in response to early Pliocene ocean warming. *Paleoceanography*, *30*(6), 789–801. https://doi.org/10.1002/2014pa002704

Haug, G. H., Ganopolski, A., Sigman, D. M., Rosell-Mele, A., Swann, G. E. A., Tiedemann, R., et al. (2005). North Pacific seasonality and the glaciation of North America 2.7 million years ago. *Nature*, *433*(7028), 821–825. https://doi.org/10.1038/nature03332

Haug, G. H., & Tiedemann, R. (1998). Effect of the formation of the Isthmus of Panama on Atlantic Ocean thermohaline circulation. *Nature*, *393*(6686), 673–676. https://doi.org/10.1038/31447

Hays, J. D., Imbrie, J., & Shackleton, N. J. (1976). Variations in the Earth's Orbit: Pacemaker of the ice ages. *Science*, *194*(4270), 1121–1123. https://doi.org/10.1126/science.194.4270.1121

Haywood, A. M., Tindall, J. C., Dowsett, H. J., Dolan, A. M., Foley, K. M., Hunter, S. J., et al. (2020). The Pliocene model Intercomparison project phase 2: Large-scale climate features and climate sensitivity. *Climate of the Past*, *16*(6), 2095–2123. https://doi.org/10.5194/cp-16-2095-2020

Hennissen, J. A. I., Head, M. J., De Schepper, S., & Groeneveld, J. (2014). Palynological evidence for a southward shift of the North Atlantic current at ~2.6 Ma during the intensification of Late Cenozoic Northern Hemisphere glaciation. *Paleoceanography*, *29*(6), 564–580. https://doi.org/10.1002/2013pa002543

Herbert, T. D., Ng, G., & Cleaveland Peterson, L. (2015). Evolution of Mediterranean sea surface temperatures 3.5–1.5 Ma: Regional and hemispheric influences. *Earth and Planetary Science Letters*, *409*(0), 307–318. https://doi.org/10.1016/j.epsl.2014.10.006

Herbert, T. D., Peterson, L. C., Lawrence, K. T., & Liu, Z. (2010). Tropical ocean temperatures over the past 3.5 million years. *Science*, *328*(5985), 1530–1534. https://doi.org/10.1126/science.1185435