Alina Găinuşă-Bogdan
email: alina.gainusabogdan@gmail.com

Didier Swingedouw

Pascal Yiou

Julien Cattiaux

Francis Codron

Francis Codron Amoc

AMOC as the key driver of the spread in Mid-Holocene winter temperature patterns over Europe in PMIP3 models

Keywords: PMIP3, Mid-Holocene, NAO, AMOC, inter-model spread

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Projections of the climate of Europe show a large spread for the ongoing century (IPCC Working Group 1 et al. 2013). For a given emission scenario, uncertainties are related to the internal variability, that can be large at the regional scale [START_REF] Deser | Communication of the role of natural variability in future North American climate[END_REF], and to the modeldependent sensitivity of the climate system to external forcings [START_REF] Hawkins | The potential to narrow uncertainty in regional climate predictions[END_REF]. In particular, changes of the atmospheric or oceanic circulation tend to vary a lot between models.

Future projections of the North Atlantic Oscillation (NAO), that plays a crucial role in the climate variability of Europe in winter [START_REF] Hurrell | Decadal Trends in the North Atlantic Oscillation : Regional Temperatures and Precipitation[END_REF], do not even agree on its sign [START_REF] Cattiaux | Opposite CMIP3/CMIP5 trends in the wintertime Northern Annular Mode explained by combined local sea ice and remote tropical influences[END_REF]. At multi-decadal timescales, the Atlantic Meridional Overturning Circulation (AMOC) also plays a key role for the climate of Europe through its massive northward transport of heat all along the Atlantic Ocean (Swingedouw et al. 2009[START_REF] Haarsma | Decelerating Atlantic meridional overturning circulation main cause of future west European summer atmospheric circulation changes[END_REF][START_REF] Jackson | Global and European climate impacts of a slowdown of the AMOC in a high resolution GCM[END_REF]. Its long-term projections are, similarly, very uncertain [START_REF] Weaver | Stability of the Atlantic meridional overturning circulation: A model intercomparison[END_REF], even though all the models agree that it should decrease in the coming century, especially when the Greenland ice sheet melting is accounted for [START_REF] Bakker | Fate of the Atlantic Meridional Overturning Circulation -Strong decline under continued warming and Greenland melting[END_REF]Swingedouw et al. 2015b). To improve our estimate of the response of the NAO or the AMOC to external forcing, paleoclimate data could provide key observational constraints [START_REF] Schmidt | Using palaeo-climate comparisons to constrain future projections in CMIP5[END_REF][START_REF] Harrison | Evaluation of CMIP5 palaeo-simulations to improve climate projections[END_REF].

Paleoclimate reconstructions are based on proxy records that have been collected for decades around the world. The main issue with proxies is that they are indirect records of climate variables, so that they could be subject to large uncertainties and biases. Nevertheless, except for the short instrumental era, they are the only observational estimates of the climate response to different external forcings. Even though there is no exact analogue of the future in the past, information on the climate response to changes in external forcings has a high relevance to better evaluate the sensitivity of climate models and possibly reduce the uncertainty concerning the on-going climate change [START_REF] Schmidt | Using palaeo-climate comparisons to constrain future projections in CMIP5[END_REF]), especially at the regional scale.

The mid-Holocene (MH: 6000 years Before Present) is a recent period where the climate of the northern high latitudes was considerably warmer than today in summer [START_REF] Fischer | Palaeoclimate constraints on the impact of 2 °C anthropogenic warming and beyond[END_REF], due to a different insolation forcing (Figure 1). This is an interesting climate period, when the Sahara is suspected to have been greener than today [START_REF] Claussen | The Greening of the Sahara during the Mid-Holocene: Results of an Interactive Atmosphere-Biome Model[END_REF]. This is notably why the Paleoclimate Model Intercomparison Project, phase 3 (PMIP3) initiative [START_REF] Braconnot | Evaluation of climate models using palaeoclimatic data[END_REF] proposed this period as a key snapshot to be simulated by the same climate models that produce climate projections [START_REF] Schmidt | Using palaeo-climate comparisons to constrain future projections in CMIP5[END_REF][START_REF] Harrison | Evaluation of CMIP5 palaeo-simulations to improve climate projections[END_REF]. A recent study proposed an improved reconstruction of climate over Europe for the MH based on a broad compilation of pollen data [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF]. The results showed a dipole of winter near-surface temperature anomalies, with warm anomalies in the north of Europe and cold anomalies in the south. The authors suggested this pattern could be caused by an atmospheric circulation anomaly similar to a positive phase of the North Atlantic Oscillation (NAO + ). Other studies also proposed such a change for mid-Holocene European climate in winter [START_REF] Guyard | Late-Holocene NAO and oceanic forcing on high-altitude proglacial sedimentation (Lake Bramant, Western French Alps)[END_REF][START_REF] Chabaud | Land-sea climatic variability in the eastern North Atlantic subtropical region over the last 14,200 years: Atmospheric and oceanic processes at different timescales[END_REF]. This hypothesis therefore provides an interesting testbed for climate models. [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] argued that PMIP3 simulations in general did not reproduced very well the reconstructed dipole in temperature in winter because they not show any NAO + -like changes in atmospheric circulation for the MH. This was also the case for the older PMIP2 dataset [START_REF] Gladstone | Mid-Holocene NAO: A PMIP2 model intercomparison[END_REF]. Nevertheless, It have been proposed that PMIP2 dataset still showed some skill in capturing the pattern of climate change over Europe for the MH (Brewer et al. 2007). While the new generation of PMIP3 models show larger changes in sea ice for MH conditions than PMIP2 models over the Arctic [START_REF] Berger | The sensitivity of the Arctic sea ice to orbitally induced insolation changes: A study of the mid-holocene paleoclimate modeling intercomparison project 2 and 3 simulations[END_REF], which may impact the response over Europe, a quantitative comparison of the PMIP3 models with winter reconstruction of temperature over Europe has not been made yet. Furthermore, a detailed analysis of the mechanism explaining the wintertime response in PMIP3 simulations for MH is also missing for the North Atlantic.

The large increase in summer insolation for MH conditions may have strongly affected the sea ice in polar regions, even in winter through inertia. This could then influence the continental climate, notably in winter where the heat flux release by the ocean is largest, as well as the insulation effect of sea ice [START_REF] Goosse | A mechanism of decadal variability of the sea-ice volume in the Northern Hemisphere[END_REF]. Such an impact was illustrated by [START_REF] Fischer | Evolution of the seasonal temperature cycle in a transient Holocene simulation: Orbital forcing and sea-ice[END_REF] in a transient simulation over part of the Holocene using the MPI climate model. [START_REF] Otto | Climate variability-induced uncertainty in mid-Holocene atmosphere-ocean-vegetation feedbacks[END_REF], using different simulations for the MH including a coupling to the ocean and vegetation or not, also showed that the ocean plays a key role for the amplitude of the temperature response of the northern high latitudes to the radiative forcing. The insolation during the MH was indeed very different from the pre-industrial period. The latitudinal gradient was strongly reduced during the summer season, which could have impacted the oceanic and atmospheric meridional heat transport [START_REF] Davis | Orbital forcing and role of the latitudinal insolation/temperature gradient[END_REF].

In the present study, we propose to analyze the response of climate models participating to PMIP3 under MH wintertime conditions, with a specific focus on the North Atlantic region.

We will assess the potential causes for the model spread in temperature response over Europe in winter and analyze the role played by changes in atmospheric circulation, oceanic circulation and radiative changes amplified by sea ice response. We will also compare these simulations, using classical skill scores, with the reconstruction of [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF]. This will allow us to determine the models that best reproduce the temperature dipole pattern over Europe. Our main result is that atmospheric circulation changes are not the main driver of the simulated patterns of temperature anomalies. These are instead linked to sea-ice anomalies, forced either by AMOC slowdown or summer sea-ice melt.

Experimental design

We analyze pairs of mid-Holocene and pre-industrial (PI) simulations performed with PMIP3 climate models. Only 12 models (Table 1) are retained, as we require at least 100 years of both PI and MH simulations to perform robust statistical analyses. When more than 200 model years are available, we use the last 200 years of the simulation. Throughout the analyses, we use data on the native grid for each model, except when multi-model statistics are calculated, in which case the model data is first interpolated on a regular 1°x1° grid for the temperature to be compared with [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] and a regular 2.5°x2.5° grid for the other variables.

Index Model Modeling Center Reference

A bcc-csm1- We focus our analysis on Mid-Holocene climatological anomalies with respect to the preindustrial period for wintertime (December-February averages -DJF), summertime (June-August averages -JJA) or annual means (ANM). Wherever suitable, we use a two-sided independent sample t-test to determine if and to what level these anomalies are statistically significant, considering the interannual variability in the PI and MH simulations, at each grid point. Where this is done, we shade in color the regions where the anomalies are significant at the 90% confidence level. For a more stringent test, we add contours corresponding to significant anomalies at the 99% confidence interval.

Temperature and circulation response MH winter temperature response in PMIP3 models

While the external forcing (insolation changes) is the same in the different models, their responses to MH climatic condition in winter over Europe are very different, as shown in Figure 2. This figure shows the significant winter surface atmospheric temperature (SAT) anomalies over Europe simulated with the 12 different climate models, the multi-model ensemble mean anomaly and the reconstruction by [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF]. The uncertainties in the pollen-based reconstruction can be quite large, therefore Figure 2n only shows those ΔSAT DJF,EU values that are larger in absolute value than the respective estimated uncertainties. We find a wide range of simulated SAT anomaly patterns between the 12 models, and some important differences from the mean reconstructed pattern. We quantify these model-reconstruction differences in Figure 3, showing the pattern correlation, pattern root-mean-square difference (RMSE; pattern RMSE obtained by removing the respective spatial average from each field before the RMSE calculation) and full-field RMSE for each model compared to the [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] reconstruction, as well as for the multi-model ensemble mean. [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] proposed that the MH temperature anomaly pattern over Europe could be caused by a NAO-like atmospheric circulation anomaly. To explore this possibility in the PMIP3 database, Figure 4 shows the statistically significant sea-level pressure differences (ΔSLP) between the climatological wintertime (DJF) MH and PI simulations over the North Atlantic (NA) sector for the 12 analyzed models. All models show significant MH winter SLP anomalies, of the order of 0.5-1.5 hPa, with corresponding 1000-hPa wind speed anomalies (ΔU horiz,DJF,1000 , not shown) of the order of 0.2-0.7 m/s. Figure 5 shows the multi-model mean circulation anomalies at different levels in the troposphere. Significant anomalies are found at all atmospheric levels. The most notable feature is a negative zonal wind speed anomaly in the subtropics (around 30 o N), maximum at upper levels (Δu 200 ) and over the western half of the North Atlantic. This feature is present in all models, but with structures varying in strength, extent and zonal orientation. We can trace it to a local weakening of the jet streams (in some models, accompanied by a slight northward shift of the eddy-driven jet) and, for most models, a more marked separation between the eddydriven and subtropical jet streams over the North Atlantic (illustrated in Figure 5a for the multimodel ensemble mean). This results in an anticyclonic horizontal wind speed anomaly and local high pressure in the West Atlantic northward of the negative Δu 200 region. For most models, this anticyclonic anomaly dominates the pattern of mean atmospheric circulation change over the North Atlantic. In those cases, this anomaly structure is found throughout the atmospheric column (as illustrated in Figure 5c for the 500 hPa level), down to the surface (Figure 5e), explaining the positive SLP anomalies found in most models over the southwestern North Atlantic (Figure 4, Figure 5d). Even the MRI model, which has no significant anomalies in the lower atmosphere, shows a similar significant anomaly in the jet stream.

MH winter atmospheric circulation anomalies in the PMIP3 models

While this positive pressure anomaly in the Western Atlantic has an equivalent barotropic structure that extends throughout the troposphere, the positive SLP anomaly seen over North Africa in the multi-model mean (Figure 5d) is on the contrary very shallow, disappearing or even changing sign at 500 hPa. This baroclinic structure is linked to a cold temperature anomaly in the atmospheric column, which can be explained by the negative DJF radiative forcing. The reduction of the radiative forcing is strongest at low latitudes (Figure 1), and has a larger impact over continental regions like North Africa.

Atmospheric circulation vs. European temperature anomaly patterns

Figure 4 shows that several PMIP3 models simulate a positive SLP anomaly over the southern part of the North Atlantic and a negative SLP anomaly over the northern part of the basin during MH winters. This pattern could be roughly construed to be NAO + -like, as [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] argued could cause their reconstructed DJF ΔSAT EU pattern. This section examines whether such a relationship between a NAO + -like circulation anomaly and a ΔSAT EU pattern can be found within our multi-model ensemble.

A key question is whether the differences between the simulated temperature anomaly patterns between models can be attributed to the differences in their simulated atmospheric circulation anomalies over the North Atlantic. More specifically, we want to assess whether the models that simulate a positive north-south temperature gradient over Europe are only the models that also simulate a NAO + -like circulation anomaly over the North Atlantic.

A qualitative comparison of Figure 2 and Figure 4 indicates that this is not the case. While the hypothesis of a NAO + -like cause of marked north-south gradients in ΔSAT DJF, EU could hold for a few models (e.g., CNRM-CM5), the ensemble of models does not support it. Some models exhibit a NAO + -like circulation without a strong north-south ΔSAT DJF, EU gradient (e.g., IPSL-CM5A-LR, MIROC-ESM), while other models show a marked ΔSAT DJF, EU gradient without simulating a NAO + -like circulation anomaly (e.g., GISS-E2-R). For some model simulations showing a negative SLP anomaly around Iceland, which can look like part of a NAO + pattern, just like over Africa the pressure anomaly disappears with height (not shown), and is most likely a thermal structure caused by the local warming at the surface.

For a quantitative assessment, Figure 6 evaluates the relationship between first-order indices of the large-scale temperature and circulation changes, namely the north-south ΔSAT DJF, EU gradient and the NAO + -like ΔSLP NA . The north-south ΔSAT DJF, EU gradient is defined as the mean temperature difference between a northern and a southern region over Europe chosen to approximately reflect the Arctic Oscillation signature on SAT EU , as shown in Mauri et al. ( 2014-top left panel in Figure 6). The NAO + -like ΔSLP NA is defined as the mean sea-level pressure difference between regions chosen to approximately reflect the centers of action of the North Atlantic Oscillation index (cf. Hurrell & Deser 2010 -bottom right panel in Figure 6). The scatterplot between the two indices in the top right panel of Figure 6 shows a lack of significant correlation that confirms that the PMIP3 models do not support a straightforward link between MH temperature anomaly patterns over Europe and a NAO + -like circulation anomaly over the North Atlantic.

Factors determining the simulated ΔSAT DJF,EU meridional gradients

Figure 2 shows that the strongest ΔSAT DJF, EU signals are found, for most models, close to the latitude boundaries of the region considered. Enlarging the domain of our analysis places the European temperature anomalies patterns in the larger-scale context (Figure 7a). We find that the modeled north-south ΔSAT DJF, EU contrasts are dominated, to the north, by warm temperature anomalies over the Nordic Seas and, to the south, by cold temperature anomalies over northern Africa.

To a large extent, the negative SAT anomalies over northern Africa (and the positive SLP anomalies) are related to the negative radiative forcing in winter at these latitudes (Figure 1). This anomaly is strongest at low latitudes and has an amplified effect on SAT over continents vs. oceanic regions, possibly due to the lower heat capacity of land than water.

Inter-model differences of ΔSAT over the southern region in Figure 5 are indeed small (0.3 o C inter-model standard deviation) relative to the ones over the northern region (0.8 o C), consistent with a lower ΔSAT inter-model standard deviations over northern Africa than over the Nordic Seas (Figure 7b). Thus, for the rest of the paper we will focus mainly on what determines the inter-model spread in ΔSAT DJF over the Nordic Seas, that dominates the north-south contrast. The spatial patterns of both the ΔSAT DJF signal in the Nordic Seas and its inter-model spread (Figure 7a,b) are strongly reminiscent of those of the DJF sea-ice cover anomalies (Figure 7c,d). Sea-ice cover changes have a strong impact on near-surface temperatures via air-sea fluxes and albedo feedbacks [START_REF] Swingedouw | Sensitivity of the Atlantic Meridional Overturning Circulation to the melting from northern glaciers in climate change experiments[END_REF](Swingedouw et al. , Årthun et al. 2017)). This influence extends to neighbouring land regions: the inter-model differences in ΔSAT DJF in land areas of northern Europe (Figure 7e) are well correlated (r 2 =0.68, p<0.01) with the differences in ΔSIC DJF over the Nordic Seas. Furthermore, these differences in ΔSAT DJF are also well correlated with the changes in latent heat flux from the ocean to the atmosphere over the Nordic Seas. Thus, we argue that it is the difference in latent heat flux released in winter by the ocean to the atmosphere and associated with the insulation effect of the presence (or not) of sea ice that explains the response of ΔSAT DJF over northern Europe land. A similar mechanism was indeed highlighted by Årthun et al. (2017) to explain the impact of Nordic Seas SST variability on the temperature of Scandinavia.

Thus, to understand the differences between the wintertime MH European temperature anomalies between the models, we mainly need to understand the origins of the differences between the simulated sea ice cover anomalies in the Nordic Seas region.

Origins of winter sea-ice decline

Atmospheric circulation anomalies

While NAO-like anomalies do not directly determine the inter-model differences in winter temperature anomalies over Europe in PMIP3 ensemble, atmospheric circulation changes between MH and PI winters could potentially control the sea ice anomalies in the Nordic Seas, through anomalous atmospheric heat transport.

To test this hypothesis, we tried to calculate the total atmospheric heat transport into the Nordic Seas region. This computation should be based on changes of both the mean and transient circulation. Unfortunately, the PMIP3 database does not routinely include data at sufficiently high frequencies to allow for the calculation of the latter term, nor does it provide estimates of atmospheric heat transport computed at the model time step. We can thus only offer a partial answer, related to the mean wind anomalies only.

Whether we consider the near-surface horizontal temperature advection due solely to the change in the mean winter near-surface horizontal wind velocities between PI and MH (estimated as , where stands for MH-PI differences, for the horizontal wind velocity and

-𝛥𝑈 • ∇𝑇 𝑃𝐼 𝛥 𝑈 𝑇 𝑃𝐼
for the temperature for pre-industrial simulation; Figure 8) or the change between the mean circulation-induced horizontal temperature advection in PI vs. MH winters (estimated as -∆

; not shown), we find no correlation of these respective anomalies with the winter sea

(𝑈 • ∇𝑇)
ice cover anomaly over the Nordic Seas region (Figure 8). More generally, two other elements of the energy budget argue against an increased poleward atmospheric transport in MH winter:

at the top of the atmosphere, the latitudinal gradient of insolation decreases (with more cooling at lower latitudes), and the ocean heat transport tends to increase (due to AMOC strengthening as shown below). Both would typically lead to a compensating equatorward energy transport by the atmosphere, away from the sea ice region. Such a response has actually been demonstrated in MH simulations run with the IPSL-CM5A-LR model [START_REF] Saint-Lu | The role of El Ni??o in the global energy redistribution: a case study in the mid-Holocene[END_REF].

We thus find no evidence that the changes in the atmospheric heat transport may be responsible for ΔSIC DJF in the models.

Summer sea ice anomalies

Another potential source of the inter-model differences in winter sea ice anomalies is the summer (JJA) sea ice response, with a memory effect between the two seasons. Indeed, the insolation is far larger in the Mid-Holocene compared to the pre-industrial period during summer (Figure 1), potentially leading to large melting of sea ice and heat storage in the ocean in summer, which may then persist into winter due to the large inertia of the ocean. There is a fine interplay between winter and summer sea ice: winter sea ice thickness heavily influences the ice pack vulnerability to summer melt, and thus has an important effect on summer sea ice cover [START_REF] Berger | The sensitivity of the Arctic sea ice to orbitally induced insolation changes: A study of the mid-holocene paleoclimate modeling intercomparison project 2 and 3 simulations[END_REF]. Conversely, summer melt allows for ocean surface heat uptake and renders new sea ice formation more difficult during autumn and winter. Thus, we can also expect a ΔSIC JJA effect on ΔSIC DJF [START_REF] Berger | The sensitivity of the Arctic sea ice to orbitally induced insolation changes: A study of the mid-holocene paleoclimate modeling intercomparison project 2 and 3 simulations[END_REF]). This latter effect may play an important role in determining the inter-model ΔSIC DJF differences. Indeed, the different models may exhibit different sea ice responses to the summer radiative forcing (Figure 1) for example through different strength of feedbacks, like the local cloud response (cf. Knudsen et al. 2015;Abe et al. 2015). These differences in summer may then project on the winter responses.

We thus check if the relationship between ΔSIC JJA and ΔSIC DJF is consistent between the models. Figure 9 indeed shows a strong multi-model relationship between ΔSIC JJA and ΔSIC DJF , both globally over the Northern Hemisphere high latitudes (r 2 =0.74, p<0.01), and regionally over the Nordic Seas (r 2 =0.73, p<0.01). While it is likely that a memory of the seasonally forced ΔSIC JJA can influence ΔSIC DJF , it is also possible that both summer and winter sea ice anomalies are influenced by a large-scale, common factor such as a change in the Atlantic meridional overturning circulation (AMOC), which we investigate in the following section.

Change in AMOC

Figure 10 shows the changes in the AMOC between the simulated MH and PI climates. All 12 models analyzed show significant anomalies in the mean AMOC, with considerable inter-model variability in terms of both anomalies (Figure 10, colors) and background state (PI AMOC -Figure 10, contours). In general, the AMOC is strengthened in MH compared to PI. This is consistent with [START_REF] Born | Sea ice induced changes in ocean circulation during the Eemian[END_REF], who proposed that for a climate with warmer summer (the Eemian in their case), the decrease in sea ice transport through Fram Strait may change the salinity in the Labrador Sea and enhance local convective activity, feeding a stronger AMOC.

AMOC MH anomalies vary between models not just in terms of intensity, but also in terms of structure, complicating their quantitative comparison. Nevertheless, maxima of the meridional stream function have been shown to be strongly related to the meridional ocean heat transport (e.g. [START_REF] Boning | Deep water formation and meridional overturning in a high-resolution model ofthe North Atlantic[END_REF]. Thus, they represent a reasonable parameter to consider in the investigation of sea ice changes. Figure 11 shows the relationship between winter sea-ice cover anomalies and anomalies in the maximum of the climatological Atlantic meridional stream function at three different latitudes:

26 o N and 45 o N, typically used when describing the AMOC strength (e.g., [START_REF] Huck | Low-frequency variations of the large-scale ocean circulation and heat transport in the North Atlantic from 1955-1998 in situ temperature and salinity data[END_REF][START_REF] Mccarthy | Measuring the Atlantic Meridional Overturning Circulation at 26 degrees N[END_REF], and 60 o N, closer to the Nordic Seas region. The DJF SIC anomalies are in all cases negatively correlated with the AMOC strength anomalies. The correlation coefficient and the slope of the corresponding linear regressions increase with latitude, so that the strongest link of ΔSIC DJF,NS with the meridional ocean circulation is found for ΔAMOC at 60 o N. The correlation coefficient is 0.89 (p-value < 0.01) and the linear regression line has a near-zero intercept, indicating that, as would be expected, an intensification of the AMOC at high latitudes is linearly associated to negative SIC anomalies in the Nordic Seas and viceversa. The correlation remains high even at 26°N (r=0.75, p<0.01), indicating that this relationship is not a purely local effect.

Changes in the subpolar gyre may also strongly affect the heat transport to the high latitudes [START_REF] Born | Multiple Equilibria as a Possible Mechanism for Decadal Variability in the North Atlantic Ocean[END_REF]. However, the barotropic stream function was not available for all the models, limiting the analysis of this effect. Since the strength of the subpolar gyre is closely related to the strength of the convection in the subpolar gyre [START_REF] Born | Multiple Equilibria as a Possible Mechanism for Decadal Variability in the North Atlantic Ocean[END_REF], we suggest that the changes in AMOC would also be reflected in subpolar gyre strength, both leading to an increase of meridional heat transport.

The changes in oceanic heat transport act on longer than seasonal time scales and may thus affect sea-ice cover in both winter and summer, likely contributing to the correlation between winter and summer anomalies highlighted in Section 6.2. Nevertheless, direct changes in summer sea ice in response to changes in insolation are also due to regional sensitivity and feedbacks that may differ among the models (e.g., [START_REF] Massonnet | Constraining projections of summer Arctic sea ice[END_REF], and so could remain a relatively independent factor from AMOC strengthening to explain the inter-model spread.

Combination of multiple factors

In the previous subsections we have analyzed one-to-one relationships between ΔSIC DJF,NS and a suite of atmospheric and oceanic variables. However, the inter-model differences in the winter sea ice anomalies are most likely explained by a combination of different factors which affect the sea ice cover in different measures depending on the model.

We have investigated a wide array of variable combinations to best explain the inter-model variations in ΔSIC DJF,NS through multiple linear regressions. We have considered parameters physically linked with the winter sea ice anomalies, like the temperature advection due to the mean circulation change in the Nordic Seas region, the surface downward radiative flux anomalies, the background state (PI) sea ice cover, the summer sea ice cover anomalies and the change in the AMOC strength.

The low number of climate simulations (12), the reductive framework (lack of sensitivity experiments) and the intrinsic high level of connectivity between climate variables through nonlinear relationships (i.e., variables not independent) limit the potential of such a multivariate analysis. Nevertheless, some information can be distilled on the basis of this exercise.

While in most combinations and for most predictor variables, the high p values associated to their respective coefficients did not allow for meaningful regressions, the p values associated to the coefficient for the AMOC strength anomaly at 60 o N were always very small (on the order of 10 -4 -10 -5 ) and the coefficients themselves were always found to be negative, indicating a robust effect of ΔAMOC on ΔSIC DJF,NS in the direction already described in Section 6.4.

One multiple linear regression did yield significant results:

with , , and an adjusted r 2 of 0.96. 𝑝 𝛼 0 = 0.04 𝑝 𝛼 1 = 10 -4 𝑝 𝛼 2 = 4 • 10 -5

Figure 12 shows the results of the multiple linear regression model in Eq. 1. Figure 12a shows the good correspondence between the reconstructed and original ΔSIC DJF,NS and Figure 12b shows, for each model, the percentage of sea ice cover change associated with each term in the regression model. According to the model of Eq. ( 1), the summer sea ice anomalies and AMOC strength anomalies have comparable effects on ΔSIC DJF,NS , with the AMOC playing a particularly important role in CNRM-CM5 and GISS-E2-R. The adjusted r 2 can be mildly improved (to a value of 0.97) by adding the pre-industrial sea ice cover over the Nordic Seas region, SIC DJF,PI, NS , as an extra predictor. However, the associated p value is 0.015, indicating a 1.5% chance that this predictor has in fact no effect on ΔSIC DJF,NS .

Furthermore, the coefficient for this extra term is positive in this regression model (indicating that climate models where more winter sea ice is present in the Nordic Seas region in preindustrial simulations tend to have positive or less negative Mid-Holocene winter sea ice anomalies) but is negative in all other variable combinations analyzed. The role, if any, of the background sea ice cover on the SIC DJF, NS anomalies thus remains uncertain.

The main caveat of the model in Eq (1) is that the two predictor variables are not independent.

ΔAMOC max,60°N and ΔSIC JJA,NS are correlated with an r 2 of 0.34, significant at the 95% level. It is indeed very likely that the AMOC anomaly affects winter sea ice both directly, through oceanic heat transport in winter, and indirectly, by influencing sea ice anomalies during summer which, in turn, through ocean-sea ice-atmosphere feedbacks, can have a lasting effect on the winter anomalies. Nevertheless, the AMOC may be affecting the annual sea ice anomalies and, since sea ice seasonality is model-dependent, adding the summer sea ice as a second predictor improves the results of the regression analysis. Such a regression could thus allow to explain some of the spread that can be explained by local sea ice feedback and not only the AMOC changes.

Summer sea ice cover anomalies are actually not entirely determined by the AMOC differences.

As shown in Figure 9, the Mid-Holocene experiences an important positive summer radiative forcing in the high latitudes and the sensitivity of each model to this forcing is very likely to have a strong effect on the summer sea ice anomalies. Since Mid-Holocene forcing is seasonal and unevenly distributed over the globe, classical definitions of climate sensitivity do not apply and could not be used to directly evaluate the model sensitivity effect on ΔSIC JJA,NS and thus on ΔSIC DJF,NS . We have attempted to use different first-order proxies of this effect that would be less correlated with the AMOC change (such as summer SAT or surface downward shortwave radiation anomalies over large domains in the Northern Hemisphere) as a replacement for ΔSIC JJA,NS in the regression model above, but with unsatisfactory results in terms of p values and correlation coefficients. the basin, in thermal balance with positive temperature anomalies. For some models these SLP anomalies could be interpreted as a NAO + signal, although the vertical structure of the NAO is much more barotropic.

While some models exhibit a NAO + -like atmospheric circulation anomaly near the surface in winter, the ensemble of 12 PMIP3 models does not support a link between this type of largescale atmospheric response and the MH temperature anomaly patterns over Europe. Instead, the simulated MH European winter temperature anomaly patterns are, in most models, determined by the radiatively forced cold anomalies at low latitudes and the more modeldependent sea ice anomalies over the Nordic Seas, with the latter responsible for approximately three times more inter-model variability than the former.

The analyses carried out to understand the differences in these sea ice anomalies between the models have revealed multiple footprints of the ΔSIC feedbacks onto other fields often dominating the climatological anomalies in the Nordic Seas region. We have shown here the signature of these feedbacks on the near-surface air temperature, but similar signatures have also been found (not shown) on other variables. A very tight collocation between strong negative (positive) sea ice anomalies and negative (positive) near-surface horizontal wind divergence anomalies was found in all models for which the latter could be calculated, likely explained by a response of the near-surface atmospheric circulation to the sea ice anomalies and to the associated air-sea flux anomalies. A local increase in cloud cover over negative ΔSIC regions (and vice-versa) has been found in most models, engendering marked local anomalies of the downward longwave and shortwave radiative fluxes at the surface (positive and negative respectively over regions of sea ice loss and vice-versa). A strong ΔSIC feedback has also been found on the surface sensible heat flux (SH). In 3 out 12 models this corresponds to a simple positive response of SH to a positive ΔSIC or vice-versa. In most models, however, the response is more complex, taking the form of a strong positive-negative SH dipole over the ΔSIC region depending on the background sea ice cover, which shows a subtle interplay of sea surface temperature and air temperature anomalies in these regions.

These ubiquitous sea ice feedback signatures confirm the importance of a correct sea ice representation in the models for representing the regional climate over Europe. All but three models show a reduced winter sea ice cover in the Nordic Seas region in Mid-Holocene compared to the pre-industrial climate. While not allowing for a quantitative model evaluation of the sea-ice cover anomaly, a recent compilation of multi-proxy reconstructions of Mid-Holocene ΔSIC which expands on the data set used in [START_REF] Klein | Model -data comparison and data assimilation of mid-Holocene Arctic sea ice concentration[END_REF] (Seidenkrantz, 2015, personal communication) qualitatively confirms this result. We can thus conclude that the models simulating an increase in NS winter sea ice cover in the Mid-Holocene, namely FGOALS-g2, MIROC-ESM and CCSM4, are unrealistic for the Mid-Holocene climate in the North Atlantic sector.

Sea ice feedback signatures dominating the climatological winter anomalies hinder the quest for the causes of SIC anomalies in the model results. Nevertheless, we were able to highlight a robust multi-model relationship between the MH winter SIC anomalies in the Nordic Seas region and the anomalies in the Atlantic meridional overturning circulation, which points to the AMOC as an important driver of the MH response of European winter temperatures.

The MH-PI differences in the atmospheric heat transport due to transient eddies form an important potential source of information on the drivers of the SIC anomalies. However, these differences could not be evaluated due to the lack of availability of the corresponding fields in the PMIP3 data base. This type of drawback will be solved thanks to the inclusion of such fields in the outputs from future model intercomparison exercises such as the upcoming CMIP6 and PMIP4.

representing the strong air-sea interaction as recently highlighted by the so-called signal to noise paradox by [START_REF] Scaife | A signal-to-noise paradox in climate science[END_REF].

The cause for the changes in the AMOC has not been analyzed in the present paper since it necessitates a dedicated study. As noticed before, this could be related to a decrease in sea ice transport through Fram straight due to a lower sea ice production when insolation is larger in summer [START_REF] Born | Sea ice induced changes in ocean circulation during the Eemian[END_REF]. Another interesting factor could be related to the influence of the atmospheric circulation on the AMOC. For instance, variability of the AMOC over the last 60 years has been often related to NAO variations [START_REF] Delworth | The North Atlantic Oscillation as a driver of rapid climate change in the Northern Hemisphere[END_REF]Swingedouw et al. 2015a). The mechanism is related to a cooling of the surface Labrador Sea during NAO + phases, which increases the Labrador deep water formation and a few years later the AMOC [START_REF] Eden | A 20-year coupled ocean-sea iceatmosphere variability mode in the North Atlantic in an AOGCM[END_REF]. Thus, if the atmospheric circulation during the MH changed towards a NAO + -like configuration, then it likely amplified the AMOC, as highlighted in Shi and Lohman (2016) for the ECHAM-FESOM-MPIOM climate model. In our present PMIP3 ensemble, the NAO changes are not the sole driver, since the relationship between our indices of NAO and AMOC changes are not significantly correlated (not shown).

If this interpretation is correct, it appears that the MH European winter temperatures can be both a testbed to estimate possible changes of the NAO and the AMOC in the near future. Since, as we find here, the inter-model differences in MH temperature anomalies are mainly explained by AMOC changes, we conclude that the MH is indeed an interesting case study to evaluate the sensitivity of the simulated AMOC to changes in radiative forcing. Nevertheless, to produce a correct assessment of models, we need to have a finer, quantitative AMOC reconstruction for this period. The use of a methodology combining climate models and available observational datasets, as employed for the deglaciation by [START_REF] Ritz | Estimated strength of the Atlantic overturning circulation during the last deglaciation[END_REF] for instance, is an interesting approach to be developed. This may help to better decipher the sensitivity of the AMOC to different forcings, which may provide new leads to reduce the very large uncertainty [START_REF] Weaver | Stability of the Atlantic meridional overturning circulation: A model intercomparison[END_REF] concerning the future of the AMOC.

Figure 1 :

 1 Figure 1: Zonally-averaged climatological monthly differences in the incoming shortwave radiation at the top of the atmosphere between mid-Holocene and pre-industrial simulations.

Figure 2 :

 2 Figure 2: (a)-(l) Significant differences (color: at 90% confidence level; inside contours: at 99% confidence level) between mean wintertime near-surface air temperature (SAT) in mid-Holocene and pre-industrial simulations run with 12 different climate models; (m) Multi-model ensemble mean difference between MH and PI climatological wintertime SAT; (n) Mauri et al. (2014) reconstruction of the MH climatological wintertime SAT anomaly -only data where the signal is larger than the estimated uncertainty is considered.

Figure 3 :

 3 Figure3: Statistics comparing the simulated and reconstructed[START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF]) mid-Holocene wintertime SAT anomalies over Europe. Correlation (values shown on the right-side vertical axis): spatial correlation coefficient between the two fields compared; Pattern RMSE: root-mean-square error when the respective spatial average has been removed from each of the two fields compared; Full-field RMSE: total root-mean-square-difference between the two fields. The model data, including the multi-model ensemble mean, are sorted from the highest to the lowest spatial correlation with the reconstructed data. The statistics are only calculated using those grid points where the reconstructed signal is stronger than the estimated uncertainty (see Figure4n).

Figure 4 :

 4 Figure 4: Significant differences (color: at 90% confidence level; inside contours: at 99% confidence level) between mean wintertime sea-level pressure in mid-Holocene and pre-industrial simulations run with 12 different climate models.

Figure 5 :

 5 Figure 5: Multi-model ensemble mean wintertime: (a) 200 hPa horizontal wind speed in preindustrial (PI), mid-Holocene (MH) simulations and the corresponding MH-PI difference; (b) 500 hPa geopotential height MH anomaly; (c) 500 hPa horizontal wind speed anomaly; (d) sea-level pressure anomaly; (e) 1000 hPa horizontal wind speed anomaly.

Figure 6 :

 6 Figure 6: Top left: PMIP3 ensemble mean anomalies of SAT (in °C) in DJF for MH minus PI. The blue box defines the northern region and the red box, the southern region used for calculating the north-shout SAT gradient. Bottom right: PMIP3 ensemble mean anomalies of SLP (in hPa) in DJF for MH minus PI. The blue box defines

Figure 7 :

 7 Figure 7: MH winter anomalies of (a, b) near-surface air temperature and (c-d) sea ice cover over the North Atlantic sector. (a) and (c) show the multi-model ensemble mean anomalies and (b) and (d) show the ensemble standard deviation. (e) Relationship between simulated near-surface air temperature over land in northern Europe and sea ice cover anomalies averaged over the region highlighted in (a) and (c) for the 12 models analyzed. (f) Same as (e) but with sea ice cover anomalies replaced by latent heat flux anomalies from the ocean to the

Figure 8 :

 8 Figure 8: Multi-model relationship between wintertime sea ice cover anomalies and mean 1000hPa horizontal circulation anomaly-related temperature advection in the NS region. The red lines indicate the linear regression and the 95% confidence interval for this regression.

Figure 9 :

 9 Figure 9: Relationship between MH winter and summer sea ice cover anomalies (a) averaged over the Northern Hemisphere high latitudes (above 60 o N); (b) averaged over the Nordic Seas (rectangular region contoured in Figure 7a, c). The text boxes in the bottom-left corner of the figures show the associated squared correlation coefficients. The red lines indicate the linear regressions and the 95% confidence interval for the regressions.

Figure 10 :

 10 Figure 10: Color: significant (at 90% level) MH anomalies of the climatological annual mean (ANM) Atlantic meridional stream function; contours: PI climatologies (contour spacing: 6 Sv; solid/dashed contours represent positive/negative values) for the 12 models analyzed.

Figure 11 :

 11 Figure 11: Relationship between mean winter MH sea ice cover anomalies over the Nordic Seas region and the MH anomaly in the maximum climatological Atlantic meridional stream function at (a) 60 o N, (b) 45 o N and (c) 26 o N. The squared correlation coefficients of the associated linear regressions are shown in the respective panels.

Figure12.

  Figure12. Results of the multiple linear regression of ΔSIC DJF,NS onto the MH anomaly of the AMOC maximum at 60 o N and ΔSIC JJA,NS : (a) correspondence of the reconstructed and original ΔSIC DJF,NS ; (b) winter sea ice cover anomaly associated to each model predictor, for the 12 models analyzed.

Table 1 : List of PMIP3 models included in this study.

 1

		1	Beijing Climate Center (BCC)	Xin et al. (2013)
			Centre National de Recherches Météorologiques/
	B	CNRM-CM5	Centre Européen de Recherche et Formation Avancées en Calcul Scientifique	Voldoire et al. (2013)
			(CNRM/CERFACS)
			Commonwealth Scientific and Industrial
	C	CSIRO-Mk3-6-0	Research Organization/ Queensland Climate	Rotstayn et al. (2012)
			Change Centre of Excellence (CSIRO-QCCCE)
	D	IPSL-CM5A-LR Institut Pierre-Simon Laplace (IPSL)	Dufresne et al. (Dufresne et al. 2013)
			State Key Laboratory of Numerical Modeling for
	E	FGOALS-g2	Atmospheric Sciences and Geophysical Fluid Dynamics/Institute of Atmospheric Physics	Li et al. (2013)
			(LASG/IAP)
	F	FGOALS-s2	LASG/IAP	Bao et al. (2013)
			Atmosphere and Ocean Research Institute
	G	MIROC-ESM	(The University of Tokyo), National Institute for Environmental Studies, and Japan Agency for	Watanabe et al. (2011)
			Marine-Earth Science and Technology
	H	HadGEM2-ES	Met Office Hadley Center (MOHC)	Collins et al. (2011)
	I	MPI-ESM	Max Planck Institute for Meteorology (MPI-M) Jungclaus et al. (2013)
	J	MRI-CGCM3	Meteorological Research Institute (MRI)	Yukimoto et al. (2012)
	K	GISS-E2-R	NASA Goddard Institute for Space Studies	Schmidt et al. (2006)
	L	CCSM4	National Center for Atmospheric Research (NCAR)	Gent et al. (2011)

  ∆𝑺𝑰𝑪 𝑫𝑱𝑭,𝑵𝑺 = 𝜶 𝟎 + 𝜶 𝟏 ∆𝑨𝑴𝑶𝑪 𝒎𝒂𝒙.𝟔𝟎 𝒐 𝑵 + 𝜶 𝟐 ∆𝑺𝑰𝑪 𝑱𝑱𝑨,𝑵𝑺

Acknowledgements

We wish to thank Marrit-Solveig Seidenkrantz, Hugues Goosse and Frédéric Hourdin for fruitful discussions that have benefitted our study. This project would not have been possible without the availability of PMIP3 data in a centralized, consistent format, for which we are deeply thankful to the climate modeling centers listed in Table 1, to the groups behind the ESGF portal, and in particular to Guillaume Levavasseur and Jean-Yves Peterschmitt, as well as the FileFinderAR5 developed by Patrick Brockmann. To analyze the PMIP3 data, this study benefited from the IPSL Prodiguer-Ciclad facility supported by CNRS, UPMC, Labex L-IPSL which is funded by the ANR (Grant #ANR-10-LABX-0018) and by the European FP7 IS-ENES2 project (Grant #312979). AGB was supported by the Swedish Research Council Grant No. C0629701 (MILEX). The research leading to this study also received funding from the French National Research Agency HAMOC project (Grant ANR-13-BS06-0003).

Discussion and conclusions

We have analyzed Mid-Holocene winter anomalies over the North Atlantic sector in 12 PMIP3 models with the purpose of comparing the simulated European winter temperature anomalies with reconstructed data, and shedding light on their drivers.

Comparing the simulated winter near-surface air temperature anomalies with the recent pollenbased reconstruction of [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] for the Mid-Holocene, we found some qualitative agreement with the paleo-reconstruction for a few models, with pattern correlation coefficients above 0.6 for 5 out of the 12 models, mostly reflecting a north-south temperature anomaly gradient over Europe. All models underestimate the amplitude of the winter anomalies, however, as well as the amplitude of this meridional gradient, compared to the [START_REF] Mauri | The influence of atmospheric circulation on the mid-Holocene climate of Europe: A data-model comparison[END_REF] reconstruction. Nevertheless, the uncertainties of this reconstruction do allow for a weaker gradient as found in the models, so that there is no striking inconsistency between this reconstruction and the models response to the Mid-Holocene forcing.

We have first analyzed the simulated large-scale atmospheric circulation anomalies over the North Atlantic sector. Our analysis has revealed a consistent climatological winter atmospheric response in Mid-Holocene consisting of two elements common to all models, though varying in spatial extent and intensity. The first is a dynamically forced anti-cyclonic anomaly over the south-western part of the North Atlantic, with a strong barotropic component. This structure is associated with a local weakening of the jet stream and a more pronounced separation between the eddy-driven and subtropical jet streams over the North Atlantic sector. The second element is a positive sea-level pressure anomaly over northern Africa with a shallow baroclinic structure. This SLP anomaly is a result of the colder atmospheric temperature in response to the strong negative radiative forcing at low latitudes, the effects of which are amplified over continental regions. In the same way, a shallow low-pressure anomaly is found in the north of Another potential key player neglected here is the feedback from vegetation changes in the high latitudes, which could also amplify the pattern found in the MH simulations [START_REF] Wohlfahrt | Evaluation of coupled ocean-atmosphere simulations of the mid-Holocene using palaeovegetation data from the northern hemisphere extratropics[END_REF], potentially bringing the model simulations closer to MH observations. We also know that PMIP3 models systematically underestimate the monsoon anomalies in the Mid-Holocene (the "Green Sahara" event, [START_REF] Claussen | The Greening of the Sahara during the Mid-Holocene: Results of an Interactive Atmosphere-Biome Model[END_REF]. It is plausible that this event was linked to a northward shift of the Hadley Cell [START_REF] Claussen | The Greening of the Sahara during the Mid-Holocene: Results of an Interactive Atmosphere-Biome Model[END_REF], which would therefore not be captured or at best be severely underestimated by the models. We conjecture that all models underestimate large-scale anomalies in the atmospheric circulation that may have taken place [START_REF] Gaetani | Understanding the mechanisms behind the northward extension of the West African Monsoon during the Mid-Holocene[END_REF] and may have affected temperature patterns over Europe during MH winters. Our study, however, highlights the potential relevance of other effects, which deserve further study.

While previous, paleo-reconstruction-based studies, often attribute Mid-Holocene winter European temperature anomalies to a NAO + -like atmospheric circulation anomaly, our analysis of PMIP3 models suggests the possibility of an additional or even alternative explanation coming from a larger geographical context. In this new view, based on the most skillful PMIP3 models, European winter temperature anomaly patterns are determined, to the south, by strong low-latitude temperature responses to the seasonal negative radiative forcing and to the north by feedbacks of the negative sea ice anomalies over the Nordic Seas on the atmosphere. Model simulations also suggest that a major driving force for these sea ice anomalies is an intensification of the Atlantic meridional overturning circulation in the Mid-Holocene compared to the modern period, as indirectly suspected in [START_REF] Masson | Mid-Holocene climate in Europe: What can we infer from PMIP model-data comparisons?[END_REF]. This is in agreement with reconstructions based on sortable silt in the North Atlantic, a proxy of deep water circulation, which indicate a stronger deep current during the MH [START_REF] Kissel | Variations in the strength of the North Atlantic bottom water during Holocene[END_REF].

The amplitude of the north-south gradient over Europe is weaker in the models than in the reconstruction. We can suggest that this is related to a general bias of climate models in