Abstract. Contributing around half of the oxygen produced on Earth, marine photosynthetic production is one of the main mechanisms for carbon fixation, with a central role in the oxygen cycle. The triple isotopic composition of atmospheric oxygen (17Δ), measured in ice cores, provides a global integrator of past biospheric oxygen fluxes, and by extension carbon fluxes. However, deconvolving the signal of 17Δ requires to isolate the oceanic biosphere productivity (17Δocean). Here, we present the first implementation of 17Δocean in the intermediate-complexity climate model iLOVECLIM. The three main processes controlling 17Δocean, i.e. photosynthesis, respiration, and air-sea gas exchange, are explicitly represented and evaluated under preindustrial conditions. Model results show overall good agreement with available measurements, particularly in the Pacific Ocean. In contrast, systematic overestimation is found in the Southern Ocean. At fixed stations, seasonality is reproduced but with underestimated amplitude. These discrepancies mainly reflect challenges in representing remineralization and oxygen minimum zones, and highlight opportunities to refine the representation of primary productivity and vertical mixing. Overall, this new implementation provides the first coupled model framework for simulating 17Δocean, both as a diagnostic of biogeochemical processes and as a tool for reconstructing past changes in marine productivity. Extending the implementation to the terrestrial biosphere will further allow reconstruction of the past global biosphere and direct comparison with 17Δ records from ice cores.

Abstract Reconstructing global temperatures from sparse proxy data remains a key uncertainty in paleoclimate science. Here, we apply an emergent constraint framework that combines statistical relationships from paleoclimate model simulations with proxy‐based estimates to infer unknown climate variables. We compiled a new multiproxy data set of tropical land temperature reconstructions for the Last Glacial Maximum (LGM), which implies tropical terrestrial cooling of −4.2°C (−5.0 to −3.5°C, 95% Confidence Interval, CI). Applying our emergent constraint framework to this terrestrial temperature reconstruction, we derived updated estimates of tropical sea surface temperature (−2.6°C; −3.3 to −2.0°C, 95% CI) and global mean surface temperature (−5.4°C; −6.5 to −4.4°C, 95% CI). These observation‐based estimates support LGM cooling toward the colder end of published ranges, implying that very low equilibrium climate sensitivity (ECS < ∼3°C) is unlikely. Our method is adaptable for different time periods and climate variables, providing a flexible framework for integrating proxy and model information.

The use of species distribution models (SDMs) in paleoecology enhances its contribution to ecology and conservation biology, such as estimating climate change impacts and informing conservation efforts. Acipenser sturio and A. oxyrinchus are closely related sturgeon species currently arranged in two disjunct distributions across the Atlantic Ocean. However, paleogenetic findings highlighted that both species were present in European waters over at least five millennia. Their joint European history during the Holocene period might have been driven by the thermal hypothesis, with A. oxyrinchus being a cooler-water species than A. sturio. Given their current poor state of conservation, quantifying their environmental preferences at a time when the two species were more abundant could help shape future conservation/restoration measures. Paleo-species distribution models were built for each species based on paleoclimate variables and archeozoological remains, most of which were genetically identified at the species level. Response curves for air temperature and precipitation were then applied to estimate relative continental paleohabitat suitability (HSI) in watersheds across four regions (i.e., the Baltic Sea, the Northern Seas, the South European Atlantic shelf, and the Mediterranean Sea) throughout the Holocene. The estimated niches confirmed the thermal sensitivity difference between the two species, with A. oxyrinchus having a narrower thermal range towards colder annual temperatures. Comparing HSI suitability between the two species, the four regions were globally all suitable for A. sturio but not for A. oxyrinchus, which showed much greater variability. Catchments with outlets to the Mediterranean Sea and the South European Atlantic shelf were mainly suitable for A. sturio. In contrast, those with an outlet in the Baltic Sea were slightly more suitable for A. oxyrinchus than for A. sturio. Catchments with outlets to the Northern Seas showed comparable suitability for both species. This superposition in near-past suitable habitats implied that the two species can potentially exist under the same environmental conditions and as such suggested future risks of biotic interactions which should be further acknowledged in potential stocking practices.

Abrupt climate changes of the last glacial period, including Dansgaard-Oeschger cycles and Heinrich Stadials (HSs), were contemporaneous with profound disruptions of the Atlantic Meridional Overturning Circulation (AMOC). However, their impacts on the western North Atlantic and adjacent North American continent remain poorly resolved. Here we present a high-resolution multiproxy reconstruction from core MD99-2203 (Cape Hatteras, 35 degrees N), integrating mineralogical and reworked dinocyst data, pollen and foraminiferal assemblages, and iLOVECLIM model simulations to identify HSs and investigate vegetation-ocean linkages during HSs 5-3. In contrast to the western European margin, each HS is marked by continental cold conditions coeval with warm sea surface temperatures (SST) and the onset of iceberg-rafted debris (IRD) deposition near Cape Hatteras, yet the regional landscape remained dominated by open forests. A close up of HS 4 reveals a distinct two-phase structure: an early stage of synchronous cooling in land and SST, associated with peak IRD deposition in the Ruddiman Belt, followed by a post-surge phase when ice-rafted debris reached 31-35 degrees N. During this phase, cold and dry continental conditions contrasted with anomalously warm sea surface temperatures, generating enhanced land-sea thermal gradients. This contrast increased humidity and, together with rising atmospheric CO2 concentration, may have triggered transient forest recovery, which coincided with subsequent surface cooling. Comparisons with transient iLOVECLIM simulations successfully capture the early phase and reveal insightful spatially heterogeneous feedbacks in the later period. Our results reveal complex, evolving land-ocean interactions during abrupt climate events, providing new insights into regional climate dynamics and impacts of AMOC weakening.

East Antarctic ice cores registered variations in local surface temperature and atmospheric composition over the past 800,000 years. However, identifying the precise timing between changes in temperature inferred from the ice and greenhouse gases records is challenging. It assumes that the modern δD–temperature relationship holds over time and there are uncertainties in the age offset between ice and air at a same depth. An alternative approach uses the δ15N of N2, measured in the air-phase, as a proxy for past temperature and snowfall accumulation. Still, a complete understanding of the processes governing δ15N changes is missing. Here we present a δ15N record from the EPICA Dome C ice core covering 800–100 ka. Comparison with a firn densification model driven by δD-based temperature and accumulation reconstructions reveals major discrepancies during glacial intervals at high obliquity. Sensitivity tests with an intermediate complexity climate model reveal that this mismatch is linked to δD not capturing the full variability of the temperature and accumulation. Our finding supports the use of δ15N as a climate tracer to make the link with greenhouse gases concentration and suggests combining δD and δ15N to best reconstruct temperature and accumulation.

Internal variations of climate can significantly influence global warming trends, especially at the continental scale, and could contribute to the recent abnormal observed warming over Europe. Model-based studies highlight that centennial variability of the North Atlantic can strongly affect this sector. However, a lack of high-resolution paleoclimate data does not allow a proper evaluation of the real existence of such a variability mode nor its amplitude. Here, we compile a series of annual proxy-based reconstructions over Europe from diverse sources to demonstrate and confirm the presence of such multi-centennial climate variability mode and quantify its amplitude. We show that this mode is closely tied to the internal variability of the Atlantic overturning circulation (AMOC) both in proxy-based reconstructions and climate models. When combined with instrumental observations, we show that the phase of this mode is crucial to be known. Indeed, results indicate that an internally-generated strengthening of the AMOC can explain a large part of the warming in the early 20th century and the relative cooling in the second half of this last century. A change in phase of this mode since the early 2000s is able to explain the observed amplified warming over Europe, which is projected to persist until the 2050s. According to an observational-constraint approach, this mode of variability could amplify the forced projected warming in Northern Europe by more than 58% in the next three decades. These results underscore the importance of considering internal climate variability when assessing regional warming trends, in order to develop consistent adaptation strategies.

The late Middle Pleistocene Transition (MPT, ~ 800-670 thousand years before present, ka) was characterised by the emergence of large glacial ice-sheets associated with anomalously warm mid North Atlantic sea surface temperatures (SST) enhancing moisture production. Still, the moisture transport across Eurasia towards high northern latitudes is poorly constrained despite its potential role as feedback mechanisms feeding the ice caps. To reconstruct late MPT moisture production and spreading, we combine records of upper ocean temperature and pollen-based Mediterranean forest cover, a tracer of westerlies and precipitation, from a subtropical drill-core collected off SW Iberia Margin, with records of East Asia summer monsoon (EASM) strength and West Pacific surface temperatures, and compare them with the iLOVECLIM model simulations. We observe that the strongest Mediterranean forest development occurred during Marine Isotope Stage (MIS) 17, centered at 700 ka, reflecting a high amount of regional winter precipitation. In contrast, MIS 19 (~785 ka), under the influence of both similar ice volume and higher atmospheric CO2 concentration, is marked by limited forest expansion indicating lower winter precipitation in SW Europe compared to MIS 17. More interestingly, the MIS 18 glacial was more forested, reflecting stronger winter rainfall, compared to the preceding MIS 19, despite that the latter interglacial was characterised by higher insolation, sea level, atmospheric CO2 concentrations and similar warm SST. The long-term increasing trend in winter precipitation in SW Europe parallels the trend of the EASM strength that reached high levels during MIS 18. The model results show high amount of winter rainfall in SW Europe and enhanced EASM (based on the modelled East Asian δ18Ocalcite and summer precipitation) for the three MISs. Similar SW European tree fraction percentages are also modelled during MIS 18 and MIS 19, as inferred from the pollen data. In contrast to the proxy data, the simulated tree fraction is the weakest during MIS 17. The simulated winter rainfall is the highest during MIS 17, but the simulated EASM is the lowest during MIS 18. This mismatch between model and proxy reconstructions could be explained by the difficulty in quantitatively estimating the forest cover from pollen data and/or the result of a feedback process that is not well reproduced in the model such as the poor prediction of the intensity and position of the oceanic moisture source despite a robust SST simulation. Here the data show that SW European winter precipitation and EASM strength reached high levels during the MIS 18 glacial. We explained that this anomalous situation was caused by nearly-continuous moisture supply from both Pacific and Atlantic oceans and its transport to higher latitudes through the westerlies, likely fueling the accelerated expansion of northern hemisphere ice-sheets during the late MPT.

Ice cores are unique archives capturing records of past temperature (through the ice isotopic composition, e.g. δD) and past atmosphere composition over the last 800 kyr. In particular, their analysis revealed that glacial-interglacial transitions, altering the Earth's climate since the beginning of the Quaternary, are associated with significant variations in the atmospheric levels of CO2 and CH4. However, comparison of past temperatures imprinted in ice-phase and atmospheric composition records imprinted in the air-phase is difficult. Indeed, the air is trapped at a depth of 50-100 m, at the bottom of the firn, where snow transforms into ice. Therefore, at a given depth, the air is always younger than the ice and firn densification modeling is needed to estimate the age difference between the air and the ice at each level. Firn densification modeling is associated with large uncertainties when it is applied to low accumulation and low temperature drilling sites of the East Antarctic plateau.An alternative approach to reconstruct air temperature directly in the air bubbles involves analyzing the isotopic composition of N2 (δ15N). Indeed, local temperature and accumulation rate evolutions affect firn thickness and hence modulated the δ15N in air bubbles trapped at the bottom of the firn via gravitational enrichment of δ15N over large glacial-interglacial transition on the East Antarctic plateau. The observation of a robust correlation between ice core records of δ15N and δD (Dreyfus et al., 2010) confirms the strong influence of local climate on the δ15N. δ15N measurements have already been applied to determine the phasing between CO2 and temperature increases over Antarctic temperature increase associated with glacial terminations. However, this strong relationship between δ15N and δD is not necessarily valid outside of glacial terminations. Here, we address the question to what extent the δ15N can be used to infer past temperatures and to study the CO2-temperature relationship, hence circumventing age uncertainties that arise when comparing ice and gas phase measurements.We first examine the δ15N record from EPICA Dome C with respect to East Antarctic climate over the last eight glacial-interglacial cycles. We use the good agreement between δD and δ15N over Termination II as a satisfactory criterion to discern when the δ15N is a reliable proxy of past temperature. Using this criterion, we assert that the correlation between δ15N and δD is robust over the past eight terminations. Focusing on the 100-300 ka BP period, we note also three intervals characterized by a weak correlation: the glacial inceptions from MIS 7e to 7d and 7a to 6e, and the MIS 6 glacial period. To explain why δ15N and δD evolutions contrast over these periods, we connect water stable isotopes with new δ15N measurements from EDC ice core and explore various snow densification scenarios yielded by a firn model under different climate conditions at the ice sheet surface. Our study permits to identify a criterion to safely use δ15N as an indicator of the past temperature in the air bubbles of the EDC ice core to study the CO2-local temperature relationship.

Reconstructions of past sea surface temperatures (SSTs) are essential for understanding long-term climate variability, yet different proxy methods can yield divergent results. In this study, we compare Mg/Ca-derived SSTs from Globigerinoides ruber sensu stricto and Trilobatus sacculifer with clumped isotope (X47) SSTs measured on G. ruber s.s. from the same core, MD96-2048 (Indian Ocean), covering the last 1.25 million years (Ma). Using the same species and samples allows minimizing ecological and environmental biases. We find that X47-derived SSTs are systematically colder than Mg/Ca-SSTs prior to 0.4 Ma, while both proxies agree well after this point. This offset is not explained by diagenetic alteration (as assessed via SEM), nor by corrections for seawater salinity, pH, or Mg/Ca composition. The Mg/Ca-SSTs from T. sacculifer are more consistent with SSTs in the older interval, but do not fully resolve the discrepancy. We found that the apparent X47-based cooling before 0.4 Ma is not supported by seawater 518O estimates or other climate indicators. Our results suggest that X47-derived SSTs may be affected by an unknown bias in older intervals, although a combination of multiple factors explored in this study could also contribute to the observed offset.

Decades of ice core research have resulted in unique isotopic data documenting changes in Greenland temperature and ice sheet evolution over the last 130,000 years. We present here a brief overview of the ice core drilling history in Greenland, the progress in the development of deep ice cores chronologies and how Greenland ice cores can bring information on past climate variations from the Last Interglacial to the last millennium and present day, including rapid variations like Dansgaard-Oeschger events. The evolution of the Greenland ice sheet in relation with the global climate system is also investigated.

Stable water isotopes are used to infer changes in the hydrological cycle for different climate periods and various climatic archives. Following previous developments of δ18O in the coupled climate model of intermediate complexity, iLOVECLIM, we present here the implementation of the 1H2H16O and 1H217O water isotopes in the different components of this model and calculate the associated secondary markers deuterium excess (d-excess) and oxygen-17 excess (17O-excess) in the atmosphere and ocean. So far, the latter has only been modelled by the atmospheric model LMDZ4. Results of a 5000-year equilibrium simulation under preindustrial conditions are analysed and compared to observations and several isotope-enabled models for the atmosphere and ocean components. In the atmospheric component, the model correctly reproduces the first-order global distribution of the δ2H and d-excess as observed in the data (R=0.56 for δ2H and 0.36 for d-excess), even if local differences are observed. The model–data correlation is within the range of other water-isotope-enabled general circulation models. The main isotopic effects and the latitudinal gradient are properly modelled, similarly to previous water-isotope-enabled general circulation model simulations, despite a simplified atmospheric component in iLOVECLIM. One exception is observed in Antarctica where the model does not correctly estimate the water isotope composition, a consequence of the non-conservative behaviour of the advection scheme at a very low moisture content. The modelled 17O-excess presents a too-important dispersion of the values in comparison to the observations and is not correctly reproduced in the model, mainly because of the complex processes involved in the 17O-excess isotopic value. For the ocean, the model simulates an adequate isotopic ratio in comparison to the observations, except for local areas such as the surface of the Arabian Sea, a part of the Arctic and the western equatorial Indian Ocean. Data–model evaluation also presents a good match for the δ2H over the entire water column in the Atlantic Ocean, reflecting the influence of the different water masses.

Ice cores are unique archives capturing records of past temperature (through the ice isotopic composition, e.g. δD) and past atmosphere composition over the last 800 kyr. In particular, their analysis revealed that glacial-interglacial transitions, altering the Earth's climate since the beginning of the Quaternary, are associated with significant variations in the atmospheric levels of CO2 and CH4. However, comparison of past temperatures imprinted in ice-phase and atmospheric composition records imprinted in the air-phase is difficult. Indeed, the air is trapped at a depth of 50-100 m, at the bottom of the firn, where snow transforms into ice. Therefore, at a given depth, the air is always younger than the ice and firn densification modeling is needed to estimate the age difference between the air and the ice at each level. Firn densification modeling is associated with large uncertainties when it is applied to low accumulation and low temperature drilling sites of the East Antarctic plateau. An alternative approach to reconstruct air temperature directly in the air bubbles involves analyzing the isotopic composition of N2 (δ15N). Indeed, local temperature and accumulation rate evolutions affect firn thickness and hence modulated the δ15N in air bubbles trapped at the bottom of the firn via gravitational enrichment of δ15N over large glacial-interglacial transition on the East Antarctic plateau. The observation of a robust correlation between ice core records of δ15N and δD (Dreyfus et al., 2010) confirms the strong influence of local climate on the δ15N. δ15N measurements have already been applied to determine the phasing between CO2 and temperature increases over Antarctic temperature increase associated with glacial terminations. However, this strong relationship between δ15N and δD is not necessarily valid outside of glacial terminations. Here, we address the question to what extent the δ15N can be used to infer past temperatures and to study the CO2-temperature relationship, hence circumventing age uncertainties that arise when comparing ice and gas phase measurements. We first examine the δ15N record from EPICA Dome C with respect to East Antarctic climate over the last eight glacial-interglacial cycles. We use the good agreement between δD and δ15N over Termination II as a satisfactory criterion to discern when the δ15N is a reliable proxy of past temperature. Using this criterion, we assert that the correlation between δ15N and δD is robust over the past eight terminations. Focusing on the 100-300 ka BP period, we note also three intervals characterized by a weak correlation: the glacial inceptions from MIS 7e to 7d and 7a to 6e, and the MIS 6 glacial period. To explain why δ15N and δD evolutions contrast over these periods, we connect water stable isotopes with new δ15N measurements from EDC ice core and explore various snow densification scenarios yielded by a firn model under different climate conditions at the ice sheet surface. Our study permits to identify a criterion to safely use δ15N as an indicator of the past temperature in the air bubbles of the EDC ice core to study the CO2-local temperature relationship.

. The EPICA (European Project for Ice Coring in Antarctica) Dome C (EDC) ice core drilling in East 19 Antarctica reaches a depth of 3260 m. The reference EDC chronology (AICC2012) provides an age vs depth 20 relationship covering the last 800 kyr (thousands of years) with an absolute uncertainty rising up to 8,000 years at 21 the bottom of the ice core. The origins of this relatively large uncertainty are threefold: (1) the δ 18 O atm ,δO 2 /N 2 22 and total air content (TAC) records are poorly resolved and discontinuous over the last 800 kyr, (2) the three orbital 23 tools are not used simultaneously and (3) large uncertainties are associated with their orbital targets. Here, we 24 present new highly resolved δ 18 O atm ,δO 2 /N 2 and δ 15 N measurements for EDC ice core covering the last five 25 glacial - interglacial transitions as well as novel absolute 81 Kr ages. We have compiled chronological and 26 glaciological information including novel orbital age markers from new data on EDC ice core as well as accurate 27 firn modeling estimates in a Bayesian dating tool to construct the new AICC2023 chronology. The average 28 uncertainty of the ice chronology is reduced from 2,500 years to 1,800 years in AICC2023 over the last 800 kyr. 29 The new timescale diverges from AICC2012 and suggests age shifts reaching 3,800 years towards older ages over 30 Marine Isotopes Stages (MIS) 5, 11 and 19. But, the coherency between the new AICC2023 timescale and 31 independent chronologies of other archives (Italian Lacustrine succession from Sulmona Basin, Dome Fuji ice 32 core and northern Alpine speleothems) is improved by 1,000 to 2,000 years over these time intervals.

Abstract. The EPICA (European Project for Ice Coring in Antarctica) Dome C (EDC) ice core drilling in East Antarctica reaches a depth of 3260 m. The reference EDC chronology (AICC2012) provides an age vs depth relationship covering the last 800 kyr (thousands of years) with an absolute uncertainty rising up to 8,000 years at the bottom of the ice core. The origins of this relatively large uncertainty are threefold: (1) the δ18Oatm, δO2/N2 and total air content (TAC) records are poorly resolved and discontinuous over the last 800 kyr, (2) the three orbital tools are not used simultaneously and (3) large uncertainties are associated with their orbital targets. Here, we present new highly resolved δ18Oatm, δO2/N2 and δ15N measurements for EDC ice core covering the last five glacial – interglacial transitions as well as novel absolute 81Kr ages. We have compiled chronological and glaciological information including novel orbital age markers from new data on EDC ice core as well as accurate firn modeling estimates in a Bayesian dating tool to construct the new AICC2023 chronology. The average uncertainty of the ice chronology is reduced from 2,500 years to 1,800 years in AICC2023 over the last 800 kyr. The new timescale diverges from AICC2012 and suggests age shifts reaching 3,800 years towards older ages over Marine Isotopes Stages (MIS) 5, 11 and 19. But, the coherency between the new AICC2023 timescale and independent chronologies of other archives (Italian Lacustrine succession from Sulmona Basin, Dome Fuji ice core and northern Alpine speleothems) is improved by 1,000 to 2,000 years over these time intervals.

The EPICA (European Project for Ice Coring in Antarctica) Dome C (EDC) ice core drilling in East Antarctica reaches a depth of 3260 m. The reference EDC chronology, the AICC2012 (Antarctic Ice Core Chronology 2012), provides an age vs. depth relationship covering the last 800 kyr (thousands of years), with an absolute uncertainty rising up to 8000 years at the bottom of the ice core. The origins of this relatively large uncertainty are twofold: (1) the δ18Oatm, δO2/N2 and total air content (TAC) records are poorly resolved and show large gaps over the last 800 kyr, and (2) large uncertainties are associated with their orbital targets. Here, we present new highly resolved δ18Oatm, δO2/N2 and δ15N measurements for the EDC ice core covering the last five glacial–interglacial transitions; a new low-resolution TAC record over the period 440–800 ka BP (ka: 1000 years before 1950); and novel absolute 81Kr ages. We have compiled chronological and glaciological information including novel orbital age markers from new data on the EDC ice core as well as accurate firn modeling estimates in a Bayesian dating tool to construct the new AICC2023 chronology. For the first time, three orbital tools are used simultaneously. Hence, it is possible to observe that they are consistent with each other and with the other age markers over most of the last 800 kyr (70 %). This, in turn, gives us confidence in the new AICC2023 chronology. The average uncertainty in the ice chronology is reduced from 1700 to 900 years in AICC2023 over the last 800 kyr (1σ). The new timescale diverges from AICC2012 and suggests age shifts reaching 3800 years towards older ages over marine isotope stages (MISs) 5, 11 and 19. But the coherency between the new AICC2023 timescale and independent chronologies of other archives (Italian Lacustrine succession from Sulmona Basin, Dome Fuji ice core and northern Alpine speleothems) is improved by 1000 to 2000 years over these time intervals.

Climate model simulations are inherently biased. It is a notably difficult problem when dealing with climate impact assessments and model-data integration. This is especially true when looking at derived quantities such as biomes, where not only climate but also vegetation dynamics biases come into play. To overcome such difficulties, we evaluate the performance of an existing methodology to correct climate model outputs, applied here for the first time to long past climate conditions. The proposed methodology relies on the ‘Cumulative Distribution Function-transform’ (CDF-t) technique, which allows to account for climate change within the bias-correction procedure. The results are evaluated in two independent ways: (i) using forward modelling, so that model results are directly comparable to reconstructed vegetation distribution; (ii) using climatic reconstructions based on an inverse modelling approach. The modelling is performed using the intermediate complexity model iLOVECLIM in the standard global and interactively downscaled over the Europe version. The combined effects of dynamical downscaling and bias correction resulted in significantly stronger agreement between the simulated results and pollen-based biome reconstructions (BIOME6000) for the pre-industrial (0.18 versus 0.44) and mid-Holocene (MH) (0.31 versus 0.40). Higher correlation is also observed between statistically modelled global gridded potential natural distribution and modelled biomes (0.36 versus 0.41). Similarly, we find higher correlation between the reconstructed and the modelled temperatures for the MH (0.02 versus 0.21). No significant difference is found for the Last Glacial Maximum when using temperature reconstructions, due to the low number of data points available. Our findings show that the application of the CDF-t method on simulated climate variables enables us to simulate palaeoclimate and vegetation distribution in better agreement with independent reconstructions.

Abstract. Stable water isotopes are used to infer changes in the hydrological cycle for different climate periods in climatic archive and numerical climate models. Following previous developments of δ18O in the intermediate complexity model iLOVECLIM, we present here the implementation of the δ2H and δ17O water isotopes in the different components of this coupled numerical climate model, and calculate the d-excess and 17O-excess. Results of a 5,000 years equilibrium simulation under preindustrial conditions are analysed and compared to observations for the atmosphere and the ocean components. In the atmospheric component, the model correctly reproduces the first order global distribution of the δ2H and d-excess as observed in the data, even if local differences are observed. The latitudinal gradient is also correctly reproduced in our model and fits previous water isotopes enabled General Circulation Models simulations despite a simplified atmospheric component in iLOVECLIM. One exception is observed in Antarctica where the model does not correctly estimate the water isotope composition, consequence of the non-conservative behaviour of the advection scheme at very low moisture content. For the ocean, the model also simulates adequate isotopic composition in comparison to the observations, except for local areas such as in the surface Arabian Sea, a part of the Arctic and West equatorial Indian ocean. Data-model evaluation also presents a good match for the δ2H over the entire water column in the Atlantic Ocean, reflecting the influence of the different water masses. Modelled δ17O and 17O-excess are also evaluated against measurements for the two components of the model and compare to another General Circulation Model.

3 δ 18 O atm δO 2 /N 2 TAC EDC depth (m) Gas age (ka BP) Uncertainty (years) Source EDC depth (m) Ice age (ka BP) Uncertainty (years) Source EDC depth (m) Ice age (ka BP) Uncertainty (years)

The end of the Middle Pleistocene Transition (MPT, ~ 800-670 thousand years before present, ka) was characterised by the emergence of large glacial ice-sheets associated with anomalously warm North Atlantic sea surface temperatures enhancing moisture production. Still, the direction and intensity of moisture transport across Eurasia towards potential ice-sheets is poorly constrained. To reconstruct late MPT moisture production and dispersal, we combine records of upper ocean temperature and pollen-based Mediterranean forest cover, a tracer of westerlies and precipitation, from a subtropical drill-core collected off South-West Iberia, with records of East Asia summer monsoon (EASM) strength and West Pacific surface temperatures, and model simulations. Here we show that south-western European winter precipitation and EASM strength reached high levels during the Marine Isotope Stage 18 glacial. This anomalous situation was caused by nearly-continuous moisture supply from both oceans and its transport to higher latitudes through the westerlies, likely fuelling the accelerated expansion of northern hemisphere ice-sheets during the late MPT.