Advances in reconstructing the AMOC using sea surface observations of salinity

The Atlantic meridional overturning circulation (AMOC) is one of the main drivers of climate variability at decadal and longer time scales. As there are no direct multi-decadal observations of this key circulation, the reconstruction of past AMOC variations is essential. This work presents a step forward in reconstructing the AMOC using climate models and time-varying surface nudging of salinity and temperature data, for which independent multi-decadal observed series are available. A number of nudging protocols are explored in a perfect model framework to best reproduce the AMOC variability accommodating to the characteristics of SST and SSS available products. As reference SST products with sufficient space and time coverage are available, we here choose to focus on the limitations associated to SSS products with the goal of providing protocols using independent salinity products. We consider a global gridded dataset and, additionally, a coarser SSS dataset restricted to the Atlantic and with a quite low spatial resolution (order of 10 degrees vs. 2 for the model grid). We show how, using the latter, we can improve the efficiency of the nudging on the AMOC reconstruction by adding a high-resolution annual cycle to the coarse resolution SSS product as well as a spatial downscaling to account for SSS gradient. The final protocol retained for the coarse SSS data is able to reconstruct a 100-year long AMOC period (average of 10.18 Sv and a standard deviation of 1.39 Sv), with a correlation of 0.76 to the target and a RMSE of 0.99 Sv. These values can be respectively compared to 0.85 and 0.75 Sv when using the global salinity surface observations. This work provides a first step towards understanding the limitations and prospects of historical AMOC reconstructions using different sea surface salinity datasets for the surface nudging.


Introduction
The Atlantic meridional overturning circulation (AMOC), as an element of the ocean circulation, is a key component of the North Atlantic climate system (Buckley and Marshall 2016). To a greater extent the AMOC variability from decadal to longer time scales relies on climate models. The longest record of direct AMOC measurements comes from the Rapid Climate Change program (RAPID) which measures the AMOC at 26 • N. It started in April 2004 and is thus slightly longer than 15 years. This is not enough for inference about decadal variability and longer time scales. Instead, this measurement has given insights on some characteristics of the AMOC such as the large amplitude of the seasonal cycle (about 6.7 Sv) driven by the wind stress over the eastern boundary (Pérez-Hernández et al. 2015) or large interannual variability with events such as the 2009-2010 dip (Bryden et al. 2014). New systems of long-term measurements of the large-scale ocean circulation at further northern latitudes have been established since then, such as OVIDE since 1993 and OSNAP since the summer of 2014 (Mercier et al. 2015;Lozier et al. 2017). But unfortunately, as the RAPID array, they only provide short term records that do not allow robust long-term analysis. Studies therefore typically turn towards historical reconstructions of the 1 3 AMOC based on other indirect variables to understand decadal variability and longer time scales.
Due to the lack of observations, the problem of the AMOC reconstruction, or even the Atlantic climate reconstruction, is approached through complementary methods. One of them is the use of direct AMOC estimates. For example, Hughes et al. (2013) showed that using near-bottom velocity and density from moorings allowed to calculate geostrophic transport from vertical gradients of ocean bottom pressure. Alternatively, satellite data have been used to compute sea level anomaly and ocean bottom pressure to estimate velocities. Frajka-Williams (2015) combined satellite measurements, cable measurements and meridional Ekman transport from winds to explain 90% of the interannual variability of the AMOC in the period April 2004-March 2014. Past AMOC reconstructions have also used indirect observation-based approaches such as fingerprints and proxies. Proxy data infer features of the Atlantic such as temperature or even parts of the ocean circulation that can indirectly be related to the AMOC magnitude (Chylek et al. 2011;Rahmstorf et al. 2015;Chylek et al. 2012;Trenberth and Fasullo 2017;Thornalley et al. 2018;Caesar et al. 2018). Additionally, sea surface temperature (SST) or the sea level off the coast of North America have been used as fingerprints of the AMOC (Duchez et al. 2016;Little et al. 2017). Alternatively, reconstructing the climate system using surface-only air-sea data and a climate model has been previously attempted (Keenlyside et al. 2008;Zhang et al. 2010). Including such information, improved the reconstruction of ocean interior estimates such as Labrador Sea water, which is directly linked to AMOC variability. A consequence of this approach is to provide predictability associated to the initialization of the internal variability. In a perfect model framework a simple Newtonian relaxation, or nudging, of surface data proved to be successful in reconstructing the AMOC Servonnat et al. 2014;Ray et al. 2015). In particular, Servonnat et al. (2014) have shown that nudging both SST and sea surface salinity (SSS) provides the appropriate surface fluxes above the right surface density thereby allowing the formation of realistic water masses, which are later transported at depth through circulation, subduction and deep convection. Servonnat et al. (2014) modulated the strength of the SST nudging through a restoring coefficient diagnosed from the observed strength of the heat flux feedback (Frankignoul and Kestenare 2002). The corresponding restoring coefficient for SSS has less physical basis and is designed to correspond to the same time scale of two months to restore a mixed layer depth of 50 m as for SST. The use of both SST and SSS was showed by these previous studies to be necessary to avoid spurious oceanic convection when SST was nudged alone in convection or near sea-ice regions. In addition, the authors avoided to nudge in sea-ice regions to avoid perturbations of the sea-ice properties. Ortega et al. (2017), also using perfect model experiments, improved this set-up further by devising a restoring term variable in time and space. More specifically, the strength of the restoring is proportional to the mixed layer depth. This improvement enables stronger nudging in deep convection regions and leads to improved reconstructions of extreme AMOC events when compared to the use of fixed restoring values (Ortega et al. 2017). These advances in the use of surface data showed good potential to reconstruct the ocean interior. Nevertheless, we can wonder if the the direct inclusion of subsurface data in assimilation scheme might provide more straightforward improvements.
The use of surface or subsurface data to reconstruct the ocean circulation is a long-lasting debate (Dunstone and Smith 2010). Including subsurface data helps constrain density gradients, which drive the ocean circulation where geostrophy applies. However, reliable and widespread subsurface data with a sufficient sampling to guide models are not available before the ARGO era (i.e. mostly after year 2000). The use of ARGO data have been shown to improve the Mercator Ocean 0.25 • global ocean analysis and forecasting system for both temperature and salinity at different depths (Turpin et al. 2015). On the other hand, not assimilating ARGO data led to large errors. Deep-ocean observations available before ARGO lack consistent temporal resolution (ocean hydrographic sections) or large-scale spatial resolution (individual hydrographic stations). Although some attempts to use these led to improved estimates of the mean ocean circulation and low frequency variations of ocean circulation and heat transport (Ganachaud and Wunsch 2000;Huck et al. 2008), but the collected data need a strong pre-treatment to be usable for climate reconstructions. Indeed, inconsistencies in the reconstruction of the ocean structure at depth arise when comparing reanalysis products (Ray et al. 2015;Karspeck et al. 2015;Born et al. 2015). To circumvent these inconsistencies and the lack of long and reliable records in the subsurface of the ocean, we here keep the approach of using only surface data to reconstruct the historical evolution of the AMOC.
For technical and historical reasons (e.g. weather forecasting), a number of observational SST reconstructions with good temporal and spatial resolutions are available since the late 19th century [i.e. HadISST (Rayner 2003) or ERSST (Huang et al. 2017)], which are mostly derived from the International Comprehensive Ocean-Atmosphere Data Set (ICOADS) covering 1662-present (Woodruff et al. 2011). In contrast, SSS observations covering the 20th century are much less numerous (Huck et al. 2008). Reconstructions based on optimal interpolation such as EN4 (Good et al. 2013) exhibit potential limitations of the spatial thinning of the profiles considered, or require further investigation of the data quality obtained from individual profiles. Recently, Friedman et al. (2017) released a new SSS data set based on in situ observations extending from 1896 to 2013 in the Atlantic (named FDS, Friedman Data Set), later extended to 2016 in Reverdin et al. (2019). Despite limitations of the record, such as low spatial and temporal resolution, the length of the dataset is an asset to provide information on the multidecadal variability of the Atlantic as well as its potential long-term trends, which are very controversial (Caesar et al. 2018). Our goal of reconstructing the multi-decadal variability of the AMOC makes EN4 and the FDS very good candidates to extend the perfect model approach of Ortega et al. (2017) to the historical context. This would allow to have reconstructions from independent datasets. Before such an attempt can be made, a first methodological step, reported here, is to understand the implication of the limited area and coarse spatial resolution nature of the FDS dataset in a similar perfect model set up as used in Ortega et al. (2017).
This study therefore aims at defining a protocol to enable the use of available multi-decadal SST and SSS observational products over the historical period in order to reconstruct the past variability of the AMOC. The model and experimental set-up is presented in Sect. 2. In Sect. 3 we explore the different impacts of the characteristics of the FDS data resolution in the reconstruction of the AMOC. Finally, the description of the most promising protocol and its implications are summarized and discussed in Sect. 4.

Model
This study is based on the IPSL-CM5A2 climate model (Sepulchre et al. submitted). This is a technically updated version of the IPSL-CM5A-LR climate model ) with improved overall model computing performance and tuning strategy to reduce a cold bias obtained in IPSL-CM5A-LR. IPSL-CM5A-LR has been used for CMIP5 and several studies aiming at understanding the natural decadal variability in the North Atlantic. A bi-decadal mode of variability is described in Escudier et al. (2013) driven by the northward propagation of temperature anomalies interacting with the sea-ice in the Nordic Seas which influence the East Greenland Current through atmospheric interactions generating finally anomalies of the opposite sign in the Labrador Sea. These will propagate northwards restarting the cycle with an opposite phase. This mechanisms was later complemented by Ortega et al. (2015) who were able to link this 20-year cycle to the westward propagation of subsurface buoyancy anomalies in the subpolar gyre. The propagation of these anomalies intensify the zonal density gradient, with an impact on the northward transport via thermal wind balance, and thus the AMOC. This 20-year cycle was also shown to be impacted by volcanic eruptions, not only in this model but in several CMIP5 models showing similar bi-decadal variability . The IPSL-CM5A-LR model was also used for decadal prediction analysis Mignot et al. 2016), as well as in the perfect model studies focusing on surface nudging techniques to reconstruct variability in the North Atlantic (Servonnat et al. 2014;Ortega et al. 2017). Swingedouw et al. (2013Swingedouw et al. ( , 2015 and Ray et al. (2015) illustrated the capacity of surface nudging techniques to reconstruct observed historical conditions with this model.
The oceanic component of the coupled model is NEMOv3.6 [Nucleus of European Modelling of the Ocean, (Madec 2008], itself composed of the ocean model OPA9 and the sea ice component LIM2 (Fichefet and Maqueda 1997). Both use the ORCA2 tripolar grid with a resolution of 2 • on average, refined to 0.5 • in the tropics, and 31 vertical levels. These components are coupled through OASIS3 (Valcke 2013) to the atmospheric model LMDZ5A ). The atmospheric component has a resolution of 96 • × 96 • resolution (i.e. 3.75 • × 1.875 • ), and 39 vertical levels. The IPSL-CM5A2 also includes the land-surface model ORCHIDEE (Krinner et al. 2005) and the bio-geochemistry module PISCES (Aumont and Bopp 2006), integrated in NEMO.
Our protocol uses a pre-industrial control simulation of 2000 years, described in Sepulchre et al. (submitted). The oceanic mean state of the last 1200 years, selected to remove initial multidecadal adjustment, is illustrated in Fig. 1. The North Atlantic shows a surface salinity pattern that broadly matches that of observations, with a maximum around 30 • N at the centre of the subtropical gyre ( Fig. 1c and see also Sepulchre et al. submitted). The strong meridional temperature gradient between 20 • and 60 • N is also represented. As shown by the maximum depth of the mixed layer in December-January-February (DJF, Fig. 1a), the model has three main sites of convection: southern Labrador Sea, south of Iceland-Irminger Sea and Greenland Sea. Mixed layer depth (MLD) variability, expressed in terms of DJF standard deviation, shows that the latter has the largest variability. The zonally averaged mass transport, as illustrated by the AMOC stream function, is characterized by a northward flow in the upper 1000 m and a return flow below 2500 m. The maximum value of the Atlantic overturning, defined as the maximum of the streamfunction below 500 m and between 30 • and 60 • N, has a value of 10.90 ± 1.27 Sv. The stream function of the model at 26 • N and 1000 m depth has values of 9.69 ± 0.94 Sv, which is significantly lower than observations in the RAPID array McCarthy et al. (2015) that give an estimate of 16.9 ± 4.4 Sv over the period 2004-2013. 1 3

Perfect model protocol
From the 1200 years of the control pre-industrial simulation, illustrated in Fig. 2 in terms of AMOC maximum evolution, a "Target period" of 100 years is defined starting in model-year 1350. Information such as SST and SSS evolution from this period are nudged in different simulations starting from different dates of the pre-industrial simulation, named Initial Conditions 1-4 (IC). The goal here is to reproduce the climate variability of this Target  period, and in particular the AMOC, using only information on SSS and SST variability over this period. In this 1350-1449 target period, the AMOC shows several peaks well outside of the one standard deviation range, which provide a testing framework for the efficiency of the reconstruction. Four initial conditions are chosen along the period as starting points of our ensemble members of our nudged simulations (crosses in the figure). These start in year 840, 1650, 1800 and 1100 with values ranging between the percentiles 40 and 89 of the 1200 years of the control pre-industrial simulation. Although we only use the initial condition, we are interested in the decadal-tomultidecadal variability imposed by the nudging, to which the initial conditions have limited impact (Meehl et al. 2014). We considered also the 100-year trajectories in the control pre-industrial simulation that follows these initial conditions as a benchmark for the reconstructions, naming them free run periods. Finally, an extra Target period (Target 2 in brown), with low AMOC variability, is considered to further test the final reconstruction protocol retained.

Background in surface nudging
The general reconstruction strategy used here builds up on previous studies from Swingedouw et al. (2013), Servonnat et al. (2014), Ray et al. (2015), and Ortega et al. (2017). Swingedouw et al. (2013) used a simple SST anomaly nudging in an historical simulation to assimilate observed SST. They then used this simulation as initial conditions to launch hindcasts and found an added skill from this procedure as compared to historical runs. The same reconstruction simulation was analysed by Ray et al. (2015) where the effect of surface SST nudging on the deeper ocean was evaluated using available observations and ocean reanalyses over the period 1949-2005. Relatively skilful subsurface reconstruction was achieved in specific regions such as the subduction in the subtropical Atlantic and below the thermocline in the equatorial Pacific. These studies were later improved in perfect model framework using not only information from the SST, but also from the SSS of the target period. Servonnat et al. (2014) and Ortega et al. (2017) showed that this can be sufficient to reconstruct AMOC modulations, in particular large variations and peaks. Following these studies, the sea-ice covered regions are excluded from the nudging to avoid spurious oceanic convection and physical disruption of the sea ice equilibrium. Target SST and SSS information is injected in the simulation using standard Haney restoring terms in the equations of temperature and salinity (Haney 1971). This additional term can be written on the right hand side of the equations describing the time evolution of SST and SSS respectively: where 0 is the mean ocean density ( kg m −3 ), C p the specific heat of the ocean ( J K −1 kg −1 ) and h the depth of the associated surface layer (m). The dots in the equation represent the other terms of the equation: atmospheric fluxes, advection, diffusion, etc. (e.g. Madec 2008 for the way these equations appear in the NEMO code).The T and S are the restoring coefficients controlling the strength of the nudging in temperature and salinity ( W m −2 K −1 and mm day −1 respectively). Following Ortega et al. (2017), these terms vary in time and space to be proportional to the mixed layer depth. They are chosen as the maximum of a background value of −40 W m −2 K −1 and − 864 mm day −1 (named cst for gamma constant) for temperature and salinity, respectively (Servonnat et al. 2014), and the adjustment of this coefficient to the value of the mixed layer depth of each gridpoint. The formula, as in Ortega et al. (2017) follows: with cst the previous restoring coefficients, MLD the depth of the mixed layer in a gridpoint (m) and h 0 a reference depth arbitrarily fixed to 50 m (m). In other words, the restoring coefficient is enhanced only when and where the mixed layer is deeper than 50 m. The model then interpolates the monthly SST and SSS from the Target period interpolate into a daily resolution, the timescale at which the nudging of SST and SSS is performed. Servonnat et al. (2014) and Ortega et al. (2017) have paved the way to historical reconstructions of the AMOC using SST and SSS nudging with a variable restoring coefficient. However, the available SSS observations present several shortcomings. SST data, like the one from ERSST (Huang et al. 2017) or HadISST (Rayner 2003), are now provided with a space resolution of around 1-2 • over the whole last century. On the other hand, the SSS of Friedman et al. (2017) presents larger uncertainties that can severely impact our reconstructions, in particular, the spatial resolution of the gridded data. Despite this drawback, this dataset has been chosen for the present tests as it is compiled from independently validated and bias corrected in situ observations. Finally, no statistical method has been used to try to improve the data, preventing a statistical reconstruction based on no real data, which may increase the uncertainty (1)

Reconstructions set-up
(2) (Friedman et al. 2017). The FDS is provided in 32 boxes of length scale of 100-1000 km covering the Atlantic between 20 • and 70 • N (Fig. 3b). The boxes are chosen to minimize missing data in each box and to maximize the coherence of low-frequency signal across each box (Friedman et al. 2017). As our focus is on the AMOC, the good coverage of the FDS in the Atlantic basin makes it a good candidate for our reconstructions. Some regions with large spatial variability are nevertheless not covered: the Gulf of Guinea, north of the Gulf Stream or the Amazon and Orinoco freshwater plumes. Finally, the data are given only as annual anomalies over the period 1896-2016 for each of the 32 boxes. In order to evaluate the potential use of this data set for historical reconstruction, we need to better understand the limitations related to the resolution and the spatial coverage of the data set.
To asses the impact of the properties of the FDS to reconstruct AMOC variations through surface nudging, we run sets of 100-year simulations in a perfect model framework using SST and SSS from the Target period under different nudging strategies (Table 1). The first set of experiments (Glob-Full, 3 members) uses the global SSS field and the SST field at the ORCA2 resolution. This set of experiments would be equivalent to using the SSS from EN4 with a global coverage. Both SST and SSS are then nudged using the variable restoring coefficient as in Ortega et al. (2017) (Eq. 3) and this ensemble of simulations is considered here as the benchmark, i.e. the best reconstruction one can expect given the protocol based on surface oceanic data only. The Atl-Full experiment tests the effect of using SSS nudging in the Atlantic Ocean only and the FDS-Full experiment in the area covered by the FDS data only (Fig. 3a, named FDS Mask). In both cases, only the spatial extent of SSS is used for the nudging is modified, but the data are considered with the ORCA2 resolution. Finally, FDS-Ave is devised to explore the impact of the coarse SSS resolution of the binned data (Fig. 3c). These values are obtained by averaging monthly SSS Target data in each of the boxes, with one unique averaged value for all the ORCA2 model grid points within each of the FDS grid-boxes. We first assess these initial experiments before moving to the more elaborated ones shown in Table 1.

Assessment of the reconstruction
We assessed the protocols comparing the ensemble mean of the nudged experiments and the Target period. For the AMOC, its magnitude is represented as the annual local  vertical maximum of the AMOC value at 48 • N. In the IPSL-CM5A2 model, this latitude corresponds to the climatological maximum of the AMOC stream function and the northward thermal-wind . We then assessed the performance of the experiments based on two metrics of the AMOC: correlation with respect to the Target period (R, Pearson correlation coefficient) and the root mean square error, RMSE, with respect to the Target period. Both values are computed over the 100-year nudging period. As an example, for the reference experiment Glob-Full, these metrics are respectively R = 0.85 * * and RMSE = 0.75 Sv (Fig. 4). The statistical significance of the correlation coefficients is computed by taking into account the serial correlation of the time series in the estimation of the number of degrees of freedom following Bretherton et al. (1999). Significance over a 99% confidence is indicated as * * . Additional physical analyses of the most skilful configurations are subsequently presented. Finally, significance of the difference in correlations is computed through a Fisher z-transform with a 95% confidence level following Fisher (1992).

Impacts of reducing the nudging area
In this section we evaluate the efficiency of the reconstruction of different oceanic variables comparing the Glob-Full, Atl-Full and FDS-Full experiments with the the target period. In this suite of experiments, the spatial extent of target SSS field is progressively reduced, but not its resolution. Starting with the sea surface reconstruction, Fig. 5 shows the spatial correlation of SST and SSS between the nudged simulations and the Target over the full 100 years of the experiments (a-c and d-f respectively). Glob-Full, Atl-Full and FDS-Full show a temperature correlation with the Target above 0.8 almost everywhere (Fig. 5a-c). The patterns of SST correlation are quite similar among the three simulations. Polar regions show the largest differences between Glob-Full and the two configurations with reduced SSS nudging (Atl-Full and FDS-Full). The Southern Ocean shows a difference of 0.7 in averaged correlation south of 60 • S (differences significant at the 99% confidence level for double sided t test). Salinity exerts a strong impact on vertical stability at these latitudes (Swingedouw et al. 2008) and SSS nudging consequently plays a key role for SST reconstruction in this region, causing a reduced correlation of the SST in the experiments with reduced SSS nudging (Atl-Full and FDS-Full). In the case of the Arctic region, a considerable part of the correlation is non-significant, making difficult to assess the impact of reducing the SSS nudging region. Notice that nudging is not done where sea-ice is present, hence the Arctic region correlation is not directly driven by the nudging procedure.
For the SSS correlation, high correlation values correspond strongly to the extent of the mask considered for each of the experiments (Fig. 5d-f). In Glob-Full the correlation is above 0.8 in most parts of the globe except in regions where sea-ice is present, since, again, we do not nudge in the presence of sea-ice. For the other configurations there are extra-tropical regions where the mean currents at the surface propagate the signal of the limited SSS nudging. In Atl-Full, for example, the Brazilian western boundary current propagates the SSS signal further South yielding significant correlation outside of the mask limit. Instead, other longitudes towards the East, i.e. 0-40 • W, show a clear latitudinal separation around 30 • S between the regions above a 0.8 correlation and lower correlations. Also, for the FDS-Full configuration the high correlation signal propagates outside the FDS mask to large parts of the Mediterranean Sea and Nordic Seas. Significant SSS correlation is also obtained in both Atl-Full and FDS-Full in the tropical Pacific and Indian basins ( ≈ 30 • S and 30 • N, Fig. 5e, f). As already shown in Servonnat et al. (2014), it is the SST nudging which allows to reconstruct some of the variability of the atmosphere and, hence, the evaporation minus precipitation (E-P hereafter) hydrological budget close to the Equator ( ≈ 10 • S-10 • N), which drives SSS in the region (Yu 2011;Zhang et al. 2012;Zheng and Zhang 2012). This latitudinal band extends poleward by advection of the mean currents of the E-P driven SSS, explaining the structure seen in the tropical regions in Fig. 5e, f outside the Atlantic. The amount of explained variance in the Pacific and Indian basin is of the order of 67% ( R 2 ) within the equatorial band ( ≈ 10 • S-10 • N) for all experiments.
The evolution of the AMOC maximum at 48 • N for the Target and the different nudged simulations is shown in Fig. 4. The three configurations show that the reconstruction is closer to the Target (higher correlation and weaker RMSE) when a larger nudging area is used for SSS. Indeed, the reference ensemble Glob-Full has the highest correlation with the Target, with a value of R = 0.85 ⋆⋆ . The correlation of the ensemble mean decreases to R = 0.8 ⋆⋆ when nudging SSS only in the Atlantic and R = 0.74 ⋆⋆ in the FDS region. This suggests that a proper representation of surface salinity in basins outside the Atlantic does influence the AMOC reconstruction, via processes that are beyond the scope of the present work. In contrast, reducing the mask extension is associated with an increase of RMSE with respect to the Target. In particular, FDS-Full, where SSS is nudged only in the FDS region, shows an increase in the spread of the ensemble compared with Glob-Full. When comparing the two ensembles, Glob-Full and FDS-Full, we find an increase in the standard deviation of the spread in the latter ensemble (0.65 Sv) as compared to the former (0.55 Sv). This difference is significant at the 99% level according to a Fisher test. This would suggest a decrease in the efficiency of our reconstruction for FDS-Full. Nevertheless, the three AMOC correlations are all high and significant, which shows that considering SSS nudging in the FDS region only is a valid way forward to reconstruct the AMOC.
To understand the impact of surface nudging on the ocean interior, we quantify the water mass transformation [WMT, following Walin (1982) and Speer and Tziperman (1992)] in the Atlantic, north of 40 • N , i.e. in a large area including the main deep convection regions of model. This diagnostic attributes water mass transformation and associated changes in ocean dynamics to the different surface fluxes, which will help us to understand the mechanisms at play. We diagnose the impact of the diabatic processes involved in the transformation of water masses by projecting the monthly surface buoyancy forcing on isopycnals between 22 and 29 kg m −3 (Fig. 6a). For this, we have grouped density values into classes of 0.2 kg m −3 interval and we quantify separately the water mass transformation (in Sv) resulting from the surface heat and freshwater fluxes of the model atmospheric fluxes (hflux and fwflux in Fig. 6a) and fluxes resulting from the nudging term alone (hfcorr and wfcorr in Fig. 6a). The values of the transformation associated with the corrections, wfcorr and hfcorr, is included in the values of the transformations due to heat and freshwater fluxes. In the case of heat fluxes, hflux is transforming water masses at densities above 27.5 kg m −3 , about one order of magnitude higher than the heat flux correction, hfcorr. The nudging terms associated to the correction are not physical and are not reflecting any processes, but just the impact of the applied correction. Figure 6a shows the averaged WMT over the 100 years of simulations for Glob-Full, Atl-Full and FDS-Full in the region north of 40 • N . For each density class, the contribution of heat and fresh water fluxes from the atmosphere is shown as anomalies with respect to the Target while the nudging fluxes are shown in full. The transformation associated to heat fluxes (dashed) is very consistent across the nudging simulations for density classes below 27.5 kg m −3 . Above this density, which is characteristic of convection regions (Fig. 6b), WMT anomalies due to heat fluxes exhibit marked differences: while Glob-Full and Atl-Full have cumulative values between − 6.5 and − 6.1 Sv, FDS-Full cumulated WMT amounts to 0.1 Sv. An extra 6 Sv of dense water is thus generated by the FDS-Full in density classes above 27.5 kg m −3 due to heat fluxes. Freshwater fluxes anomalies contributing to WMT are comparatively smaller in all of the experiments, varying between − 0.2 and 0.2 Sv above 27.5 kg m −3 . Combining these two terms we obtain an increase WMT for FDS-Full as compared to the other experiments, mainly driven by heat fluxes in the shown density classes. Higher density transformation implies the formation of excessive dense water in density classes not present in Glob-Full nor in Atl-Full. This increased dense water formation contributes to the positive anomalies of AMOC observed in FDS-Full during some periods of the 100 years of simulation in Fig. 4.
The freshwater flux corrections arising from the nudging procedure (wfcorr, dotted) shows an increase of dense watermasses transformation from − 2 and − 3 Sv in Glob-Full and Atl-Full compared to − 0.3 Sv in FDS-Full (cumulative values above 27.5 kg m −3 ). Note that the freshwater flux correction is computed only north of 40 • N and it is associated to the SSS nudging area, which is reduced in each experiment. The reduction in the SSS nudging region could explain these differences in the freshwater flux corrections. Instead, as we have nudged the SST globally in the three experiments, the heat flux correction extension is maintained in the three experiments and the restriction of north of 40 • N to compute the WMT defines the same region. FDS-Full shows the largest correction for densities above 27.5 kg m −3 with a cumulative value of 17.7 Sv  Fig.6a). Glob-Full and Atl-Full show similar lower values of WMT around 8 Sv (peak of 1.7 Sv at 27.8 kg m −3 in Fig.6a). Variations of these corrections are associated to the variations in the region of SSS nudging. By reducing the SSS nudging region, we increase the region where we nudge SST only. Only nudging SST corrects the SST through heat fluxes and, hence, modifies the watermasses. However, since there is no SSS-nudging in FDS-Full, this transformation is unconstrained reaching higher values when compared to the ensembles where SSS nudging is active.
We located spatially the excess of water mass formation associated to heat fluxes by detecting the model gridpoints whose time mean values of fluxes and corrections in FDS-Full are larger than the 99-percentile of fluxes and corrections of the Glob-Full over the whole region. FDS-Full shows anomalously strong heat fluxes in grid points located near convective regions outside the SSS nudging region (Fig.6b). This means that only SST nudging is applied at these points. Hence an excess of dense water is formed through heat fluxes in FDS-Full in regions where SSS is not nudged. This supports the idea that although heat fluxes are promoting the creation of denser watermasses in these points, the lack of SSS nudging allows salinity differences (with respect to the target) to trigger this process. In fact, when we nudge the SSS as in Glob-Full or Atl-Full, we do not have this excessive WMT. This result is also supported by previous work from Servonnat et al. (2014) which concluded that nudging both SST and SSS is required to avoid spurious convection at high latitudes.

Impacts of the resolution of salinity data
The available SSS observations from the FDS data do not have the spatial resolution of ORCA2, the spatial resolution used in FDS-Full. Thus, in spite of its spatial restrictions, FDS-Full experiment is still an ideal case and a similar experiment cannot be produced using observed data. In order to test the impact of this limitation for future reconstructions using real data, we evaluate it still using perfect model framework. To do so, we first averaged the monthly SSS Target data in the FDS grid boxes as shown in Fig. 3c. The crude averaged data, named FDS-Ave, and the original Target, averaged over the 100-year period, differ considerably in regions of strong salinity gradients (Fig. 7a). These errors reach up to values of ±1 psu in the subpolar region. These anomalies contribute to a degradation of the reconstruction of the AMOC maximum at 48 • N (brown line in Fig. 8) with values of R = 0.31 ⋆⋆ and RMSE = 1.86 Sv . As shown in Fig. 9a, the WMT anomaly with respect to the Target, and due to heat fluxes, is severely degraded with anomalies reaching 10 Sv in most of the density classes above 27 kg m −3 .
To avoid these degradations, we construct a new synthetic SSS Target field using a convolution in which the interannual evolution of the FDS data in each of the large gridboxes is added to a seasonal cycle of SSS at the original ORCA2 resolution. The rationale is that, for historical data, although the past interannual variability of SSS is largely unknown throughout the 20th century, robust SSS climatologies are now available at relatively high resolution, at least higher than that of the FDS data [e.g. ARGO (Gaillard et al. 2016) and more recently satellite data from SMOS (Font et al. 2013)]. In practice, for the perfect model framework used here, we first create from the Target synthetic SSS monthly anomalies similar to the ones obtained in the FDS. Secondly, we compute a monthly climatology from the 100 years of the Target period ( SSS climat ) to obtain the climatology at ORCA2 resolution. Then the climatology at each individual ORCA2 grid-point is added to the interannual anomalies averaged in each grid-box. This convolution can be described mathematically as follows: for m in month, y year and (i, j) gridpoints of the ORCA2 grid in Box n of the FDS boxes.
Using this convoluted approach and repeating the same diagnostic shown in Fig. 7a, the previously obtained SSS biases are now strongly reduced (Fig. 7c). Nevertheless, absolute errors accumulated over the 100-year period remain (Fig. 7d). Such an error can be expected as we are not reproducing perfectly the variability of the SSS target.
The AMOC reconstruction at 48 • N obtained with the convoluted SSS field (Conv-Full, green curve in Fig. 8) shows an improved agreement with the Target as compared to FDS-Ave ( R = 0.61 ⋆⋆ ). The downside of this configuration is a large RMSE (reaching 2.35 Sv), which arises essentially from an overestimation of the AMOC. Although in most of the periods (specially during the high AMOC peaks of the target) we are overestimating the AMOC, the Conv-Full is also able to return to AMOC values similar to the target (around year 30, 50 and 80). The WMT analysis shown in Fig. 9a (green curve) shows an overestimation of WMT in the 27.8-28kg m −3 density class. These density ranges are similar to the ones identified in FDS-Full, which pointed out that nudging only SST in regions of deep convection could be one of the potential causes to this excess of surface water mass transformation. To avoid this situation, we further improve the protocol by limiting the SST nudging spatial extent. In the northern hemisphere, we limit SST nudging to 60 • N except in the Atlantic region where we follow the northern boundary of the FDS mask. In the southern hemisphere SST is not nudged south of 60 • S . The ensemble performed with this configuration, Conv-60NS (not shown in Fig. 8), yields an improved AMOC reconstruction with a correlation coefficient of R = 0.74 ⋆⋆ and RMSE = 1.23 Sv. This configuration, in particular, reduces the excess of deepwater formation due to heat fluxes in the 27.8-28kg m −3 density class (cyan in Fig. 9a). Heat flux corrections due to the nudging and heat fluxes of Conv-60NS are compared with the ones in Conv-Full in Fig. 9b. The marked red crosses and black dots on each grid-point show respectively an excess of heat fluxes and heat flux correction of the Conv-Full ensemble when compared to the 99-percentile of the mean heat fluxes and heat flux correction of the Conv-60NS. Figure 9 confirms that the decrease of dense water formation in the 27.8-28kg m −3 density class from Conv-Full to Conv-60NS is associated with a reduction of the excessive heat fluxes Further examination of fluxes and corrections in the Conv-60NS experiment showed that the water flux correction has a distinctive pattern associated with the structure of the boxes in the Convoluted salinity (not shown). In particular, we observed excessive corrections in the Conv-60NS at the edges of the FDS boxes when compared to the Glob-Full experiments. The potential reason is discontinuities related to the construction of the convoluted salinity, in which each FDS box is constructed individually with discontinuities at each box border. Hence, we improved further the target SSS field by smoothing it with a spatial filter. The filter is based on a weighted average of each model grid point and the nine neighbouring ones, giving an equal weight to the central grid point and the neighbouring ones. The spatial filter is used to compute a new SSS field used in the ensemble Conv-60NS-Sm. The mean water flux corrections in Conv-60NS-Sm improve the water flux correction anomalies located near the edges of the FDS boxes (not shown). With this improvement we obtained our best candidate for the reconstruction of the AMOC (Fig. 8, blue curve). The metrics of the reconstruction for the AMOC at 48 • N are indeed very similar to the best we could possibly expect given the spatial extent of the FDS data , namely FDS-Full, R = 0.76 ⋆⋆ and a RMSE = 0.99 Sv in Conv-60NS-Sm versus R = 0.74 ⋆⋆ and a RMSE = 1.0 Sv in FDS-Full.

Analysis of the final protocol
We have based the choice of our best protocol on the variability of the AMOC maximum at 48 • N . In this section we examine further the characteristics of the obtained reconstruction.

AMOC variability
We first look at the meridional coherence of the final AMOC reconstruction by considering depth-latitude AMOC averages of the subpolar (45-60 • N) and subtropical (20-45 • N) regions in the North Atlantic. We find that the subpolar region is correlated with the target to 0.70 ⋆⋆ in Conv-60NS-Sm, compared to 0.69 ⋆⋆ for FDS-Full and 0.85 ⋆⋆ for Glob-Full. In the subtropical area the correlation slightly improves for Conv-60NS-Sm with a R = 0.54 ⋆⋆ compared to the R = 0.49 ⋆⋆ for FDS-Full (both compared to 0.57 ⋆⋆ for Glob-Full). The correlations between Conv-60NS-Sm and FDS-Full are not statistically significant different (Fisher z-transform) in both the subtropical and subpolar region. Hence, Conv-60NS-Sm represents the subpolar and subtropical regions similarly to the representation expected by limiting the salinity nudging to the spatial extent of the FDS data (FDS-Full). Instead, compared to Glob-Full, we obtain a significant different degradation in the AMOC correlation in the subpolar region, not in the subtropical. This shows the importance of SSS at high latitudes for an improved representation of the subpolar overturning through nudging. It also shows that variations in the SSS extension have a limited impact in the reconstruction of the subtropical overturning. A potential improvement would be to include wind nudging following Ortega et al. (2017), as wind circulation in the Atlantic region is a key driver of the AMOC variability at subtropical latitudes (Bryden et al. 2014).
Since the variability of the Target has a strong multidecadal variability, we further test the Conv-60NS-Sm protocol that presents a different target period with lower variability: Target period 2 (years 1555-1655, brown period in Fig. 2). In Fig. 10 we show the new target period and the corresponding reconstruction (brown), together with the original Target period and the associated reconstruction with Conv-60NS-Sm (blue). For the new target, the obtained reconstruction of the AMOC has a correlation of R = 0.71 ⋆⋆ and RMSE = 0.76 Sv (Fig. 10). Differences in the correlations of the protocol with the 2 target periods are not statistically significant (Fisher z-transform). For both experiments, we also computed 21-year window correlations with respect the corresponding target to evaluate the reconstruction of decadal variability. We considered this analysis between year 11 and year 90 to maintain the length of the window. Years with significant 21-year window correlation, above 95% confidence level, are represented as dots in Fig. 10. The reconstructions provide with a range respectively of 70 and 71 significantly correlated years out of the possible 80. Conv-60NS-Sm of Target and Target 2 take respectively 10 and 20 years to have significant values. This indicates that the effectiveness of the nudging in the initial years depends not only on the initial conditions, but also on the initial variability of the target. With the same initial conditions, we have an improved efficiency of the nudging when the target has a higher variability (Target).

Regions of deep water formation
The production of deep dense water in convection regions is a key component of the overturning circulation. Here we asses the density of the water column in these regions (coloured boxes in Fig. 7c) and over the 100 years of simulations. We show density anomalies with respect to the Target averaged in the region.
In the Greenland Sea, there is a weakening of density in Conv-60NS-Sm when compared to the Target (Fig. 11a). In the upper 500 m, and still within the standard deviation of the Target, the anomaly reaches − 0.02 kg m −3 whereas at depth (500-3000 m) Conv-60NS-Sm has an anomaly outside the STD with a mean of − 0.01 kg m −3 . We hypothesis that this negative anomaly is due to the fact that the FDS mask does not reach the Greenland Sea, so that neither SST nor SSS are nudged in this region in Conv-60NS-Sm. This contributes to reduce the efficiency reproducing convection in this region. In configurations where we nudge only SST in the region, such as FDS-Full and Conv-Full, we have significantly higher density in the deep ocean (significance estimated by a two-sided t test comparing the time series of the ensembles averaged between 1000 and 3000 m). As Servonnat et al. (2014) showed, nudging SST-only can create spurious convection increasing density in the deep ocean.
In the convection region South of Iceland, there is a good match in the deep ocean between Conv-60NS-Sm and the Target (Fig. 11b). The main difference is in the upper 1000 m where the averaged negative anomaly reaches − 0.01 kg/m 3 . It originates both from differences between the convoluted SSS and the Target SSS, and Fig. 10 Reconstruction of the AMOC at 48 • N using the configuration Conv-60NS-Sm for the target with lower variability (Target 2, thick black line) and the original target (Target, thin black line). For each reconstruction, we show the mean (lines in colours) for the Target 2 (brown) and Target (blue) and and envelope limited by the standard deviation of the ensemble members at each time step. Significance of the correlations of the whole time series between each experiments ensemble mean and corresponding target are indicated on each label as * * for a 99% confidence levels. Dots indicate a 21-year running correlation significant at the 95% confidence level from an excess of surface WMT in high-density classes (Fig. 7a). All the other experiments, using directly the salinity field from the Target, namely FDS-Full and Glob-Full, remain within the standard deviation of the Target. Conv-Full has the highest density anomaly in the deep ocean outside the standard deviation. This configuration has imperfect salinity that could deliver anomalous SSS, compared to the Target, and thus modifying the superficial density. Additionally, the northern part of the convection region is only nudged to SST, promoting the formation of dense water in the wrong density classes as shown by excessive heat fluxes in Fig. 9b.
Finally, in the Labrador Sea, Conv-60NS-Sm also exhibits a lower density than the Target in the upper 500 m with an averaged value of − 0.03 kg m −3 (Fig. 11c). These anomalies decrease slightly at depth ( − 0.003 kg m −3 on average between 500 and 3000 m), reaching values much closer to the Target than the other simulations. As Conv-60NS-Sm, Glob-Full anomalies remain within the standard deviation of the Target whereas Conv-Full and FDS-Full have anomalies at depth larger than the standard deviation (FDS-Full only below 2000 m). We hypothesize that the high density at depth for Conv-Full in this convective region is developed by advection of high density at depth from both Greenland Sea and South of Iceland. In fact the main differences affecting density between Conv-Full and Conv-60NS-Sm are the limitation of SST nudging and the smoothing of the SSS field. The latter should not have major impact at convective regions to generate large differences between these two configurations. Instead, looking at the densities from the other sites, the advection of density anomalies at depth could explain the excess in density anomalies at depth in the Labrador Sea for Conv-Full.
From all these analyses, we conclude that Conv-60NS-Sm is the best approach to reconstruct the AMOC within this perfect model framework.

Conclusions
We have presented here a protocol to reconstruct the North Atlantic Ocean circulation from surface observations, with a particular focus on the AMOC. The protocol is developed in a perfect model framework using sea surface salinity and temperature nudging within the IPSL-CM5A2 model. As our final goal is to reconstruct the 20th century, we have used the spatial features of the recent SSS reconstruction from Friedman et al. (2017), FDS, and of standard SST products such as ERSST (Huang et al. 2017) or HadISST (Rayner 2003). We explore more specifically issues related to the impacts of limited coverage and spatial resolution of observed SST and SSS. The protocol is based on previous results (Ortega et al. 2017) in which the nudging of surface data uses a variable restoring coefficient which depends on the mixed layer depth in order to maintain a relaxation time-scale in the mixed layer of about 2 months everywhere (Frankignoul and Kestenare 2002).
Unlike SST, SSS observations are scarce. The FDS SSS observation-based product (Friedman et al. 2017) is one of the longest datasets over the historical period, but it is provided in 32 coarse boxes of length scale of 100-1000 km covering the Atlantic between 20 • S and 70 • N. The data set provides a regular observing record over the whole period avoiding spurious variability obtained with irregular observing periods in various reanalysis products . Using this SSS product to nudge a model presents two main issues: a reduced spatial domain and a coarse resolution. These issues have been explored here in a perfect model framework choosing a Target period of 100 years from a pre-industrial control run. The steps of our investigation are summarized in Fig. 12, which showed the correlation coefficient and RMSE of reconstructed AMOC variations at 48 • N in the various nudged simulations compared to the Target. We start from Glob-Full, a configuration where SST and SSS are nudged globally at the resolution Fig. 11 Density anomalies with respect to the Target in each of the Convective regions for the main experiments: Glob-Full (yellow), FDS-Full (red), Conv-Full (green) and Conv-60NS-Sm (blue). Standard deviation (STD) of the density in the Target simulation is shown for reference (dashed black lines). The density anomaly axis is zoomed in below 700m (note different scales) of the model (this configuration could be associated with a reconstruction using the EN4 dataset as SSS target). From this configuration, the area of SSS nudging is restricted to the FDS region (FDS-Full) while keeping the grid resolution of the Target nudging data. This does not significantly alter the performance of the AMOC reconstruction. Then, directly averaging the Target data in the coarse FDS boxes (FDS-Ave) severely degrades the AMOC reconstruction, as the fine SSS gradients within these boxes are not captured adequately. A water mass transformation analysis (WMT) shows that the details of the deep convection in the model cannot be properly captured any more. This issue is resolved by building a synthetic Target SSS field defined as a convolution between a "high-resolution" annual cycle and the coarse resolution interannual anomalies coming from the FDS resolution (Conv-Full). To avoid nudging just the SST at high latitudes, as this generates spurious convection, the extent of the nudging region for SST is finally restricted between 60 • S and at the northern limit, a mixed limit between 60 • N and the northern boundary of the FDS mask in the North Atlantic, which extends locally further north than 60 • N . Additionally, the Convoluted dataset is smoothed to reduce discontinuities in the edges of the FDS boxes leading to the Conv-60NS-Sm experiment which exhibits the same level of correlation and RMSE as the FDS-Full. This final protocol is therefore able to reconstruct the AMOC at 48 • N with imperfect salinity at the same level as having perfect salinity with the same regional limitation.
The next logical step is to apply this protocol to the actual SSS data from Friedman et al. (2017) in the historical context, with the goal of reconstructing the AMOC from the early 20th century. As the FDS only provides annual means and the current nudging protocol was explored using monthly data, a time convolution method will also have to be used to account for seasonal non-linearities. Both FDS and the SSS from EN4 will be considered to perform these reconstructions. These two datasets also have an associated uncertainty which has not been explored in this work. This uncertainty, together with the provided SST uncertainty in standard datasets, could could influence the performance of the protocols to reconstruct the AMOC. Hence, an assessment of the uncertainty will be a crucial step required to ensure the quality of our reconstructions ). These issues are currently explored and will be reported along with the first historical reconstructions in a followup study. In this future work we will be able to compare the obtained results with other ocean reanalysis products presented in Carton et al. (2019) and previous AMOC reconstructions and estimates of the AMOC in the 20th century (e.g. Balmaseda et al. 2007;Caesar et al. 2018;Thornalley et al. 2018;Rahmstorf et al. 2015). Finally, this work will also be used as initial conditions to the decadal prediction of the IPSL-CM6-LR model. We expect that the improved initialization of the oceanic component, by the protocol presented here, will provide a better decadal prediction skill (Counillon et al. 2014;Polkova et al. 2019).
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.

Fig. 12
AMOC correlation and RMSE with respect to the Target for each experiment. The arrows indicate the different steps taken from one nudging protocol to the next to reach to final protocol (Conv-60NS-Sm)