Resolving the controls of water vapour isotopes in the Atlantic sector

Stable water isotopes are employed as hydrological tracers to quantify the diverse implications of atmospheric moisture for climate. They are widely used as proxies for studying past climate changes, e.g., in isotope records from ice cores and speleothems. Here, we present a new isotopic dataset of both near-surface vapour and ocean surface water from the North Pole to Antarctica, continuously measured from a research vessel throughout the Atlantic and Arctic Oceans during a period of two years. Our observations contribute to a better understanding and modelling of water isotopic composition. The observations reveal that the vapour deuterium excess within the atmospheric boundary layer is not modulated by wind speed, contrary to the commonly used theory, but controlled by relative humidity and sea surface temperature only. In sea ice covered regions, the sublimation of deposited snow on sea ice is a key process controlling the local water vapour isotopic composition.


S
O undergo isotopic fractionation during phase transitions of water. Therefore, they can be used as integrated tracers of hydrological processes in the atmosphere. Their relative abundances compared with H 16 2 O, expressed as δ 18 O and δ 2 H (see the Methods section), have been measured and used for many applications in climate-related studies, e.g., as proxies for past temperature 1 and precipitation 2,3 changes, variations of atmospheric moisture source conditions and transport pathways 4,5 .
During phase changes, equilibrium and kinetic fractionation processes differently affect δ 18 O and δ 2 H. The deuterium excess 6 , hereafter d-excess, has been defined to quantify the kinetic effects (see the Methods section), such as those occurring during oceanic evaporation 6 or snow formation from supersaturated vapour at low atmospheric temperatures 7 . Merlivat and Jouzel 8 , hereafter referred to as MJ79, developed a first theoretical model of isotope fractionation processes during evaporation from the ocean surface, which is still widely used. Applying their theoretical concept to the Earth's global water cycle, MJ79 introduced the so-called "closure assumption", assuming an equality of the isotopic composition of the net evaporated flux and the initial moist air above the ocean surface. According to this model, the strength of the d-excess signal in vapour is related to the relative humidity of the near-surface air with respect to the saturation vapour pressure at the ocean surface (RH sea ), as well as to the sea surface temperature (SST). The theoretical considerations by MJ79 led to different interpretations of past d-excess variations recorded in polar ice cores. They have been used as proxies of changes of the moisture source relative humidity 9 or SST [10][11][12] . The latter interpretation requires the assumption of negligible relative humidity variations during past climate changes, which has been recently challenged 13 . Regionally limited water vapour isotopic observations document a primary influence of the relative humidity on dexcess variability, while the influence of SST remains difficult to assess in this context [14][15][16] . Based on the evaporation theory and observations at the microphysical scale at the atmosphere-ocean interface, the model by MJ79 also considers an impact of wind speed on kinetic fractionation processes during evaporation and subsequently on the d-excess in the atmospheric vapour 8,[17][18][19] . The importance of this wind-speed effect could not be validated so far by vapour d-excess observations performed only in coastal stations like Bermuda or Iceland 14,15 , and is therefore still under debate.
In polar regions, variations in sea ice extent are supposed to affect regional precipitation amounts 20 and to be reflected in the water isotopic composition [21][22][23] . Current understanding of the impact of sea ice on the vapour isotopic composition is, however, still limited by the number of observations available in sea ice covered areas 24,25 .
Here, we present a unique new dataset of ship-based in situ isotopic measurements of vapour and ocean surface water, conducted with an identical instrumental setup over 2 years for a large range of oceanic surface conditions at the basin scale of the Atlantic and Arctic Oceans, contrary to previous measurements more confined in area and time [25][26][27] . Our measurements, together with theoretical calculations and atmospheric simulations, allow for the first time the assessment of the variability of the water isotopic signal on the first order (the δ values of isotopic abundancies for different species) and on the second order (the dexcess signal) under various climate conditions. For the process of oceanic evaporation, our dataset is consistent with the role of meteorological RH sea and SST, but rules out the theoretically assumed influence of wind speed on the d-excess of the initial vapour. Furthermore, the sublimation of snow deposited on top of sea ice is identified as a crucial process determining the nearsurface vapour isotopic composition in sea ice-covered areas.

Results
Spatial and temporal variations. Our observations, recorded onboard of the research vessel Polarstern, cover the period 29 June 2015 to 30 June 2017 and extend over a large range of latitudes in the Atlantic sector (i.e., Atlantic Ocean and the Atlantic regions of the Arctic and Southern Oceans), from the North Pole in the Arctic to the Weddell Sea in coastal Antarctica (see the Methods section for details).
All atmospheric measurements and simulated values are presented at a 6-h temporal resolution (see Fig. 1 Outputs from an atmospheric general circulation model with an explicit diagnostic of stable water isotopes (a so-called isoGCM) nudged to meteorology (see the Methods section) are compared with the observations. The atmospheric isotopic measurements are very well reproduced by one of the simulations (ECHAM final ) for the complete observational period (see Fig. 1). The data-model agreement is excellent for temperature (correlation coefficient R 2 = 0.96, Pearson correlation p-value < 0.01, for N = 2389 data points), specific humidity (R 2 = 0.97, p < 0.01, N = 2513), δ 18 O (R 2 = 0.82, p < 0.01, N = 2466) and δ 2 H (R 2 = 0.85, p < 0.01, N = 2466) and good for the second-order d-excess signal (R 2 = 0.46, p < 0.01, N = 2466).
Deuterium excess controls during oceanic evaporation. First analyses focus on the data obtained over open-ocean regions without land or sea ice upwind, within the observation period 29 June 2015 to 1 July 2017 (see the Methods section for details on the selection criteria). The corresponding dataset (N = 1070 simultaneous vapour isotopic and meteorological observations) is distributed along the ice-free Atlantic region, from 81°N, near Svalbard, to 74°S in the Amundsen Sea (see Supplementary  Fig. 3). Thus, the isotopic dataset and the related climate variables, e.g., SST and RH, are subject to both spatial and temporal (synoptic, seasonal and interannual) variations. Theoretical calculations derived from the MJ79 model and outputs from the isoGCM are used here to evaluate the evaporation processes influencing our observations. A large range of meteorological conditions is covered by the observations from the open ocean, with very dry (+ 41%) to supersaturated (+ 125%) RH sea values and SST from −1.8°C to + 29.1°C. The correlation between both parameters is very low (R 2 = 0.10, p < 0.01). In our observations, the d-excess values in near-surface vapour are anti-correlated with RH sea and correlated with SST (with R 2 = 0.62 for RH sea and R 2 = 0.50 for SST, p < 0.01). A multivariable linear regression of d-excess against both parameters indicates that constant d-excess values are distributed along oblique lines in an RH sea /SST diagram (see Fig. 3). We obtain the empirically estimated function (with R 2 = 0.76 and a root-mean-square error of 3.4‰) with d-excess in ‰, RH sea in %, and SST in°C. We note that this empirically estimated function includes both spatial and seasonal temporal changes in the evaporation conditions. While it covers a large range of meteorological conditions in the Atlantic sector, caution should be taken to apply this empirical function to isotopic data sampled in very different spatial or temporal domains.
To investigate the influence of the wind speed on the d-excess signal, we first focus on the distribution of observed d-excess values against wind speed ( Supplementary Fig. 4). To filter out the primary control of RH sea and SST on the d-excess signal, we sort our observational dataset into several categories, where both RH sea and SST vary within a small range only. In each of these categories, no relationship can be observed between the wind speed and the d-excess values. Under the assumption that the measured d-excess values are caused by kinetic fractionation occurring during local oceanic evaporation, our results indicate that these fractionation processes are independent of the concurrent wind speed.
Next, we compare our open-ocean water isotopic measurements with calculations of the atmospheric boundary layer water vapour isotopic composition, based on the MJ79 evaporation model. In this model, an influence of wind speed on the kinetic fractionation during evaporation is considered, as wind will affect the surface roughness by generating waves, which in turn might alter the evaporation flux. Based on laboratory experiments 19 , the MJ79 model assumes a smooth and a rough wind regime (below and above 7 m s −1 surface wind speed, respectively) with distinct kinetic fractionation coefficients for both regimes. Three different parameterisations of the kinetic fractionation are applied in our calculations. In the first parameterisation, a discontinuity is assumed in the kinetic fractionation coefficients at the windspeed threshold of 7 m s −1 , as suggested by MJ79. The two other parameterisations use constant kinetic fractionation coefficients, identical to those applied either below or above this wind-speed threshold (see details in the Methods section). Our calculations imply a local closure assumption by considering only the local variations of RH sea , temperature, wind speeds and oceanic surface water isotopic composition, neglecting any potential mixing of local vapour with advected air masses or convection.
The δ 18 O and δ 2 H values of all calculations are comparable, but do not match the observations (see Supplementary Fig. 5). The calculations always underestimate the short-term isotopic variations (related to synoptic variability) and overestimate the average isotopic levels compared with the observations. This overestimation can be explained by the applied local closure assumption, as the mixing of locally evaporated moisture with advected humidity is neglected in this model approach. For the MJ79 model, the closure assumption is in general not valid at the local scale, but at the global scale only 29   In contrast to δ 18 O and δ 2 H, the atmospheric boundary layer variations in the d-excess signal are primarily controlled by kinetic fractionation processes occurring during evaporation. The wind-speed-related parameterisation of kinetic fractionation in the MJ79 model strongly affects the calculated d-excess values (see Fig. 4, Supplementary Fig. 6 and 7). Different kinetic fractionation coefficients used in the M79 calculations below or above the 7 m s −1 threshold lead to different slopes in the distributions of d-excess versus RH sea for the two wind regimes (see Fig. 4). However, in our measurements, the distributions of d-excess against RH sea are nearly identical for both wind regimes (with 63% of wind-speed conditions above 7 m s −1 , see . The calculation based on the kinetic fractionation coefficients of the smooth wind regime leads to an RH sea /d-excess distribution (Fig. 4) with a similar slope as the observations, but is biased towards higher d-excess values (slope of calculated versus measured values m = 0.94, R 2 = 0.71, but with a + 4.9‰ offset, p < 0.01, see Supplementary Fig. 7). We note that none of the parameterisations are able to reproduce the lowest measured d-excess values, corresponding to the highest RH sea values (see Fig. 4). Despite the overestimation of the firstorder isotopic signals δ 18 O and δ 2 H in the local closure assumption, our observed d-excess variability can thus be reproduced by the MJ79 model approach, even on a local scale. The observations are better reproduced if constant kinetic fractionation coefficients are applied, and the best match between our data and the MJ79 model is achieved for the constant kinetic fractionation coefficients of a rough wind regime.
In the observations, the d-excess/RH sea distribution is characterised by a slope of −0.39‰ % −1 (R 2 = 0.64, p < 0.01; see Fig. 5). In the MJ79-based calculations, this distribution is slightly different from the observations for any parameterisations of the kinetic fractionation, but the deviations are smaller when using the coefficients of a rough wind regime as compared with using the ones of a smooth wind regime. For the coefficients of a rough regime, the slope is of −0.32‰ % −1 (R 2 = 0.69, p < 0.01; see Our results are based on observations within the boundary layer,~30 m above the skin layer, at which the evaporation takes place. On the opposite, the wind speed-dependent kinetic fractionation parameterisation of MJ79 is based on wind tunnel experiments performed for a limited range of wind speeds and wave types 19 . From our analyses, we cannot make any conclusive statement about the validity of the model, as we did not performed comparable (laboratory) experiments directly above the water surface. However, our results clearly indicate that the MJ79 model should not be applied in its original form for the calculation of isotopic changes in atmospheric vapour well above the ocean, e.g., as done in current isoGCMs. The wind and wavetype range investigated for the MJ79 model approach might not necessarily represent the diversity of surface oceanic conditions observed at sea. For example, a rough ocean surface with high waves might also be caused by swell, and does not have to be directly linked to high wind speeds occurring at the same time. Based on our new dataset, we rather suggest to modify the MJ79 model and use constant kinetic fractionation coefficients instead of wind-speed-dependent values. This conclusion is supported by the recently reported lack of wind-speed influence on the water Influence of sea ice on water vapour isotopic composition. Next, we focus on a subset of data gathered at high latitudes, where sea ice has been surrounding Polarstern (for local sea ice fractions higher than 0%, see the Methods section for details). The corresponding dataset contains measurements from both the Arctic and Antarctic regions (N = 854 simultaneous vapour isotopic and meteorological observations). A recent study 24 postulated an anti-correlation between vapour d-excess and local sea ice fraction, based on near-ocean-surface vapour isotopic measurements, conducted over the course of approximately 3 summer days in the western Arctic. This anticorrelation was linked to the meteorological conditions at the sea ice margin. Our measurements cover a substantially longer time period and a larger spatial scale within both the Arctic and Antarctic sectors. They do not confirm such anti-correlation but rather indicate a positive correlation, with a d-excess increase of 14‰ from open-ocean conditions to a complete sea ice coverage (see Fig. 6, Supplementary Fig. 9). However, the correlation between vapour d-excess and sea ice fraction is weak (R 2 = 0.35, p < 0.01) and does not improve when separating Arctic from Antarctic data for the analysis. A decrease of RH sea is observed with increasing sea ice coverage, related to an air temperature decrease, while the SST values cannot be lower than about To identify the potential cause of this effect, we compare our measurements to two isoGCM simulations (N = 840 comparison points with observations; see Fig. 1 and Supplementary Fig. 10). In the first simulation, we assume that the isotopic composition of a bare sea ice surface is equal to the isotopic composition of the ocean water just beneath the sea ice, which is the usual procedure in such isoGCM simulations. For the second simulation, we assume that the isotopic composition of the sea ice surface is a function of the isotopic composition of a snow layer deposited on this surface (see the Methods section for details). Sublimation from the sea ice surface to the lowest atmospheric model layer is allowed in both cases, without considering any isotopic fractionation. In the first simulation, the modelled variations of δ 18 O and d-excess are small and do not agree with the measurements (R 2 = 0.14, p < 0.01 for δ 18 O; R 2 = 0.00, p > 10 −1 for d-excess, respectively; see Supplementary Fig. 8). In the second simulation, the measured low δ 18 O and high d-excess values of vapour over sea ice-covered areas are better simulated (see Fig. 1 and Supplementary Figs. 8 and 10). Spatial and temporal variations of both parameters are reproduced (R 2 = 0.6 for δ 18 O, R 2 = 0.35 for d-excess, p < 0.01, see Supplementary Fig. 8) for measurements in both hemispheres. We conclude that the snow accumulated on top of sea ice, which has depleted δ 18 O and δ 2 H, and high d-excess values compared with seawater, is a potential additional key source determining the atmospheric boundary layer vapour isotopic composition in sea icecovered regions. We note that the applied parameterisation of the fraction of sea ice covered by deposited snow (see the Methods section) is based on a subset of our observational data and thus does not represent a strict independent proof for the importance of snow sublimation as a source for the isotopic composition of the vapour. We rate it as a first-order approach to include snow on sea ice for 15 20 Y = 0.14 X -6.82 R 2 = 0.35; p = 7.1e-82 Fig. 6 Sensitivity of the deuterium excess in vapour with respect to different sea ice coverage conditions. Relation between d-excess measurements and sea ice coverage around Polarstern for the period 29 June 2015 to 30 June 2017 (blue dots). Only d-excess measurements for sea ice coverages higher than 0% are considered. A linear regression is displayed in black. All plots are based on data with a 6-h temporal resolution future isotope modelling studies. Further observational data are certainly necessary to validate and improve this parameterisation, e.g., to take the flushing of the snow by seawater in fragmented sea ice areas into account, as well as potential isotopic fractionation effects during the sublimation of the snow.
During August 2016, measurements on the research vessel have been performed in an area with only a partial sea ice coverage in the vicinity of the Greenland ice sheet, close to the outlet of the Nioghalvfjerdsbrae glacier (latitude 79°N). Very depleted isotopic values of near-surface vapour measured during this period (δ 18 O reaching a local minimum of −37.7‰ close to the values observed at NEEM on top of the Greenland ice sheet 30 ) are matched by both isoGCM simulations, independently of the parameterisation of sublimation above sea ice. Advection of isotopically depleted vapour from the Greenland ice sheet towards the research vessel could create a signal overprinting the local vapour isotopic composition. The model simulates sublimation over Greenland with the same surface source for both parameterisations, contrary to the sublimation taking place on the sea ice, and would provide the same isotopic composition in both simulations. However, such influence of katabatic winds on our dataset is generally limited, both around Greenland and Antarctica. Within the sea ice-covered area, air masses originating from coastal regions, as filtered for the open ocean, represent~8% of the dataset. Our results concerning the sea-ice influence on d-excess do not change when filtering data points potentially influenced by such continental air masses (not shown here).

Discussion
Our results are based on direct isotopic measurements, on calculations applying the MJ79 model and on results from complex isoGCM simulations. Our measurements support the fundamental theory of kinetic fractionation by MJ79 8 concerning the influences of both relative humidity and temperature at the atmosphere-ocean interface on the atmospheric boundary layer d-excess of vapour over the oceanic surface. However, contrary to this theory, our data suggest that the kinetic fractionation is not modulated by wind speed. Considering constant fractionation coefficients with values for a rough wind regime yields best agreement between observed and modelled d-excess values. The general relationship we obtain for the distribution of d-excess as a function of relative humidity and SST is based on a compilation of observations from various climatic regions, ranging from the tropics to high latitudes. For the calculation of this relationship, we neglected the potential influence of advected moist air on our measured data. Thus, the relationship should be used with care for oceanic regions, where moisture advection might substantially contribute to the boundary layer water vapour content. For sea ice-covered regions, our results indicate that sublimation of snow on sea ice might be a key additional process, controlling the isotopic composition of the boundary layer water. This vapour can subsequently influence the isotopic signal of polar precipitations.
Hence, our results have, among others, the following implications for paleoclimate studies based on water isotope records, e.g., derived from ice cores and speleothems, as well as for present-day hydrological studies. Firstly, the variations of d-excess should be interpreted as a mixed proxy for both relative humidity and SST conditions at the moisture source, but not as a proxy for wind speed. In this regard, a 10% increase in RH sea would reduce the dexcess by~3‰, while a 10°C increase in SST would raise the dexcess by about 3‰. Secondly, at high latitudes, isotopic variations in near-surface vapour are strongly influenced by evaporated ocean water, but potentially also by a snow cover on the sea ice, which has an isotopically different source signal than ocean water. Combined with the decrease of relative humidity towards sea ice-covered areas, this leads to an~1.2‰ decrease in δ 18 O and 1.4‰ increase in d-excess for every 10% increase in sea ice coverage. This sea ice effect in δ 18 O, δ 2 H and d-excess may have an imprint on the subsequent water isotopic composition of precipitations. It may then contribute to explain, for instance, some abrupt variations of the d-excess signal recorded in Greenland ice cores at the end of the last glacial period 31 or to validate a hypothesis of past sea ice retreat at 128 ka around the West Antarctic Ice Sheet 32 . Water isotopic variations in ice cores may also be used as a proxy for regional sea ice extent in the Arctic and Antarctic sectors, in combination with other chemical proxies 33 . For this purpose, it is needed to carefully evaluate the moisture source locations of the sites where the ice cores are retrieved, as well as potential post-depositional processes affecting d-excess values in the firn layer. Another implication of our results concerns the parameterisation of future isoGCM simulations focusing on polar regions, which should explicitly consider the snow on top of sea ice, identified as a potential additional sublimation source affecting the isotopic signal. The first-order parameterisation deduced from our observational data and isoGCM experiments might be used for such modelling studies, but further independent observational data and simulation results are certainly required for improving this parameterisation.

Methods
Meteorological observations. Routinely measured meteorological data from Polarstern are used in this study. The related sensors are located at different heights: wind speeds and wind directions are measured at 39 m above sea surface, relative humidity (RH air ) and temperature (T air ) at 29 m above sea surface and water temperature (SST) at 5 m below sea surface. Air pressure (P) is measured at an altitude of 19 m, but expressed at sea level. The calibrated and validated datasets are available at a 10-min averaged temporal resolution on PANGAEA Open Access library 34 and have been averaged at a 6-h temporal resolution in this study.
The relative humidity of the near-surface air with respect to the saturation vapour pressure at the ocean surface (RH sea ) is not directly measured but has been approximated 27,35 from the observed relative humidity RH air at 29 m, corrected by the ratio of specific humidity at saturation between the temperatures at this elevation and at the sea surface (T air and SST) where q sat (T) is the specific humidity at saturation for a given temperature T and q sat (SST) is calculated for seawater at salinity 35 PSU 36 . For intercomparison with other sea surface water vapour isotopic measurement campaigns 27 , this calculation is performed using the air temperature and relative humidity corrected from the altitude at 10 m. Sea ice coverage. Due to a lack of continuous and quantitative sea ice observations during the different Polarstern cruises, the sea ice coverage surrounding the research vessel has been derived from ERA-interim reanalyses 38 at 0.75°× 0.75°s patial and 6-h temporal resolution. The sea ice coverage at a specific Polarstern location is assumed to be equal to the value of the surrounding grid cell. This dataset has been compared and is coherent with values extracted in the same manner from daily sea ice coverage data from the AMSR2 instrument on-board the GCOM-W1 satellite at a 6.25-km resolution.
Definition of deuterium excess. The deuterium excess values are computed based on the commonly used definition 6 : Water vapour isotopic composition. A Cavity Ring Down Spectroscopy (CRDS) analyser (model L2140-i, Picarro, Inc.) has been running continuously on-board of the research vessel Polarstern since the 29 June 2015, recording humidity mixing ratio, δ 18 O and δ 2 H values of water vapour at a temporal resolution of~1 s. The ambient air inlet for this instrument is located at 29 m above the sea level, connected to the analyser through an~25 -m-long tubing heated at 65°C. The humidity mixing ratio is converted into specific humidity measured by the CRDS analyser (q CRDS ) and corrected by a linear function derived from the direct comparison with specific humidity values derived from the meteorological observations (q meteo ) on-board the Polarstern, during the complete measurement period from 1-h resolution datasets: q meteo = 0.75 × q CRDS − 0.17 (R 2 = 1.0, p < 0.01, for N = 17592). q meteo is calculated based on RH air , T air and P. The precision of specific humidity measurements is estimated at 0.1 g kg −1 from the comparison of both datasets. For instrument calibration of the isotopic values, a custom-made system is used, vaporising water isotopic standards injected in liquid form and mixed with dry air provided by high-pressure gas cylinders. Four different liquid isotopic standards are used, covering the range of the expected ambient air values (δ 18 O values between −7.8‰ and −40.7‰). Recommendations for long-term calibration of CRDS water vapour isotopic analysers were followed 40,41 . Therefore, our system allows two types of calibration. Firstly, the concentration dependence 34,35 of the raw isotope measurements is corrected. Secondly, repeated corrections of the deviation of the measurements from the VSMOW-SLAP scale 42 are performed by the computation of calibration curves based on the measurements of the water standards, thereby allowing correction of the instrumental drift.
The humidity-concentration dependence of the isotope measurements is corrected based on the measured isotopic composition of each four water standards over a range of different humidity values. The results of the calibration measurements are presented in Supplementary Fig. 1. The temporal stability of this correction has been evaluated by successive measurements of this so-called humidity-response function at different times. No significant drift of this response was observed for any of the four standards, neither for successive measurements over a week (not presented separately in the graphics) nor for measurements separated by several months over the complete observational period (as shown by the different measurement sequences in the graphics). The humidityresponse function is thus considered constant in time. It does however depend on the isotopic standard used. A humidity-response function is computed for each isotopic standard as the interpolation of the distribution of all experiments with a polynomial function of fourth order. The correction of the humidity-concentration dependence for a specific near-surface vapour measurement is determined by the linear interpolation of the two humidity-response functions from the closest surrounding isotopic standards at the isotopic value of the measurement.
Calibration curves are applied to the raw data to correct deviations from the VSMOW-SLAP scale. These calibration curves are calculated based on the repeated measurement of every liquid standard for 30 min every 25 h (a standard measurement sequence consists of the successive measurement of all four calibration standards). To avoid any memory effects, averaged values and standard deviations of the standard isotopic composition are computed over the last 15 min of each injection only. Several filtration and correction steps (summarised in Supplementary Table 1) are applied to these standard measurements before computing the calibration curve. All measurements are corrected for the humidity-concentration dependence. To account for the difference in the isotopic composition of the same liquid standard stored in different bottles and used on separate injection lines, we define an arbitrary reference standard among both bottles and correct the measured isotopic values from the difference between the known isotopic value of both bottles. We remove measurements with average H 2 O values (H 2 O) below 5000 ppm or higher than 28000 ppm and standard deviations of H 2 O, δ 18 O or δ 2 H (noted σ(H 2 O), σ(δ 18 O) and σ(δ 2 H), respectively) higher than 2500 ppm, 1.5‰ and 5‰. We compute a first 14-day running average and eliminate all measurements that deviate from this running average by more than 1.5‰, 5‰ and 8‰ for δ 18 O, δ 2 H and d-excess, respectively. The observed variabilities of these selected measurements of all liquid standards are shown in Supplementary Fig. 2.
The calibration curves are calculated every time a standard measurement sequence has been performed, based on a new 14-day running average of the previously selected liquid standard measurements. These values are compared with the theoretical values of the reference standards at the time of the standard measurement sequence. If values of at least 3 standards are available, a linear regression of the measurements against the theoretical values gives the calibration curve. Otherwise, as found in the literature for such type of analysers 14 , we correct the calibration curve from the drift of the running average of the standards, which have been correctly measured and use an interpolated value of the slope between the closest calculated calibration curve.
Based on the uncertainty of both corrections from the concentration dependence and deviations from the VSMOW-SLAP scale, the measurement accuracy is estimated at 0.16‰, 0.8‰ and 2.1‰ on δ 18 O, δ 2 H and d-excess. The measurement precision on 1-h averages, estimated from the standard deviation of calibration standard measurements at a constant humidity level, is of 0.24‰, 0.7‰ and 2.7‰ on δ 18 O, δ 2 H and d-excess for humidity levels above 5 g kg −1 . It deteriorates with lower humidity levels, reaching 0.5‰, 1.9‰ and 5.9‰ for δ 18 O, δ 2 H and d-excess, for humidity levels of 1 g kg −1 . The dataset presented in this study has been averaged at a 6-h temporal resolution.
Surface water isotopic composition. Isotopic composition of the surface oceanic water has been measured from daily taken samples since 30 June 2015, filled in narrow-mouth low-density polyethylene 20-or 30-mL plastic bottles, sealed with Parafilm M and stored at +4°C from the end of the expedition until the measurement. Measurements are done with IRMS and equilibration technique at the isotope laboratory of AWI Potsdam 43 (with an accuracy better than 0.1‰ and 0.8‰ for δ 18 O and δ 2 H). For comparison with other parameters, an interpolation of this dataset has been used at a 6-h resolution.
MJ79 model under the closure assumption. For all the observations performed above the open ocean, we compute the corresponding theoretical water vapour isotopic composition in the atmospheric boundary layer over an open ocean based on the MJ79 model under the closure assumption. We assume it to be equal to the isotopic composition of the evaporation flux 8,35 . We thus express the boundary layer vapour isotopic ratio (R BL ) as a function of the surface seawater isotopic ratio (R SW ), taking both equilibrium and kinetic fractionation coefficients α eq and α k and RH sea into account: We use skin temperature at the air-sea interface from the OSTIA dataset to determine α eq values 44 . R SW is determined by the interpolated values of the isotopic composition measured in daily sampled surface oceanic water. We use three different parameterisations for the kinetic fractionation coefficients dependency on wind speed. In the first simulation (named MJ79 ref ), the kinetic fractionation coefficients (α k; 18 O and α k; 2 H for H 18 2 O and 1 H 2 H 16 O) for a smooth or a rough wind regime are used for wind speed, respectively, below or above the threshold of 7 m s −1 . We respectively apply kinetic fractionation coefficients 35,45 : α k; 18 O =1.0060, α k; 2 H = 1.0053 for a smooth wind regime; α k; 18 O = 1.0035, α k; 2 H = 1.0031 for a rough wind regime. We use the values of the measured wind speed on Polarstern (at 39 m above sea surface). In two additional simulations (respectively named MJ79 smooth and MJ79 rough ), the kinetic fractionation coefficients are set constant and independent of the wind speed, either to the value of a smooth or a rough wind regime.
Atmosphere general circulation model with water isotopes. In this study, isoGCM simulations are performed with the ECHAM5-wiso model 46 with a horizontal grid size of~1.1×1.1°(T106 spectral resolution) and 31 vertical levels. The model is nudged to ERA-interim surface pressure, temperature, vorticity and divergence fields 47 to ensure that the simulated large-scale atmospheric flow is modelled in agreement with the ECMWF reanalysis data on all analysed timescales during the years 2015-2017. For each time step of 6 h, isoGCM simulation results of near-surface vapour amount and its isotopic composition are extracted from the model grid cell encompassing the position of Polarstern. In the vertical direction, this grid cell extends from the surface to~60 m above the surface. Two different ECHAM5-wiso simulations are performed, all covering the period from January 2015 to July 2017 after a 12-month spin-up period.
In the first simulation (named ECHAM exp ), different kinetic fractionation coefficients during evaporation over open water are applied depending on wind speed: constant coefficients for a smooth wind regime, and wind speed-dependent coefficients for a rough wind regime are used for wind speeds below or above the threshold of 7 m s −1 , respectively 48 . Over sea ice-covered areas, bare ice is prescribed with an isotopic composition of ocean surface water, based on a global gridded data compilation of δ 18 O in seawater 49 .
In the second simulation (named ECHAM final ), constant fractionation coefficients for δ 18 O and δ 2 H, suggested for a rough wind regime 35 , are applied under all different meteorological conditions for evaporation processes. Over sea ice-covered areas, a 2 -cm-deep snow layer pad is assumed on top of any sea icecovered grid-cell fraction, to account for accumulation and sublimation of snow on sea ice. The isotopic composition of the bare sea ice is assumed to be equal to the composition of surface ocean waters, neglecting a potential small fractionation process occurring during the formation of sea ice. The prescribed surface ocean δ 18 O and δ 2 H values are taken from a reference global gridded dataset compilation 49 . The isotopic composition of the snow layer is determined by the isotopic composition of the accumulated snowfall. During sublimation processes, no fractionation of the snow is assumed. This treatment of snow on sea ice as a singlelayer bucket model is equivalent to treatment of snow on land surfaces in ECHAM5-wiso. The deposited snow in the model is locally controlled, without taking any advection of sea ice or snow drift into account. In reality, the isotopic composition of the sea ice surface will not only be determined by the isotope signal of bare ice or snow on top of the sea ice, but also by further processes altering the sea ice surface. For example, sea spray or breaking waves might substantially alter the isotopic composition of a snow-covered sea ice surface, especially for regions with only a minor area fraction covered by sea ice. These effects will lead to a further mixing of the isotopic signal of the original fallen snow with the isotopic composition of the surrounding ocean surface waters. To account for such processes in ECHAM5-wiso, the isotopic composition of the sea ice surface is assumed as with δ as δ 18 O or δ 2 H, δ sea ice surface the isotopic composition of the sea ice surface, δ snow bucket the isotopic composition of the snow bucket on sea ice, δ ocean the isotopic composition of the surrounding ocean surface water and f as the fraction of sea ice in each grid cell. This empirical formula is based on a comparison of measured and simulated δ 18 O and d-excess values for the period from August to October 2015 and applied for the data-model comparison over the whole measurement period of this study (see Fig. 1) afterwards.
Data filtering. For analysing evaporation processes occurring over an open ocean, without any potential influence of land-or sea ice-based processes, we filtered all the data with sea ice or land situated upwind of each measurement. The upwind area is defined by a 40°angle centred around the wind origin and a maximum distance of vapour transport within 14 h previous to the measurement, which is determined by the measured Polarstern wind speed. Only if this area is free of both sea ice and land (0% sea ice index and no land area), we consider the corresponding measurements as influenced by surface processes over the open ocean. Vice versa, the isotope data corresponding to sea ice-covered areas are selected by including all measurements, where the ERA-interim sea ice coverage in the grid cell surrounding the research vessel Polarstern is higher than 0%. The locations of all filtered datasets are displayed in Supplementary Fig. 3a, b.

Data availability
All presented instrumental and modelling data of this study are available on the PANGAEA database 50 .

Code availability
The code of the ECHAM model can be retrieved from the Max-Planck-Institut für Meteorologie and is subject to a license. The isotope enhanced version is available by personal contact to the authors.