1 Introduction

The North Atlantic is one of the regions where near-term climate predictions are most promising (Kirtman et al. 2013). Such near-term climate predictions, on interannual to decadal timescales, have a strong potential to influence our society with benefits to agriculture (Hammer et al. 2001), energy supply strategies, adaptation to global change, etc. However, these applications depend on the accuracy and reliability of the predictions (Slingo and Palmer 2011). In turn, the latter depend on a careful assessment of prediction uncertainty. Indeed, in a perfect and therefore reliable prediction system, prediction uncertainties and forecast errors are expected to be equal on average (Palmer et al. 2006). For lead times shorter than a few decades, internal variability and model imperfections have been shown to be the major contributors to the climate projection uncertainty in contrast to the uncertainty arising from emission scenarios for greenhouse gases (Hawkins and Sutton 2009). Near-term climate prediction experiments strive to reduce the projections uncertainty by carefully initialising the climate system (Meehl et al. 2013). However, even for small errors in the initial state, a large uncertainty may arise from the non-linearity of the system (Lorenz 1963). This source of uncertainty is usually taken into account by performing ensemble predictions with slightly perturbed initial conditions.

Several ensemble generation techniques based on atmospheric perturbations only, extending from random perturbations (e.g. Griffies and Bryan 1997; Persechino et al. 2013) and shifting atmospheric state by a few days (e.g. Collins and Sinha 2003; Collins et al. 2006; Yeager et al. 2012), to more elaborated methods designed to generate optimal initial perturbations, such as atmospheric singular vectors (e.g. Hazeleger et al. 2013) and breeding vectors (e.g. Ham et al. 2014), have been used for decadal predictions and predictability analyses. Although, all of these methods generate ensemble spread in the whole climate system, they neglect uncertainties in the ocean initial state that need to be taken into account at seasonal and decadal timescales. This may result in insufficiently dispersive ensembles, thereby leading to overconfident and therefore unreliable forecasts (e.g. Ho et al. 2013). Despite these generally accepted ideas, the inclusion of ocean state uncertainties in the initial ensemble spread remains challenging.

Germe et al. (2017) compared the impact of random atmospheric perturbations vs oceanic perturbations mimicking random oceanic uncertainties and found that the latter have the same impact on the future evolution of the ensemble as atmospheric-only perturbations after the first 3 months in the IPSL-CM5A-LR climate model. However, Du et al. (2012) showed that oceanic perturbations arising from different assimilation runs do affect the ensemble spread of oceanic-related variables. This latter result can be accounted for by the differences between initial oceanic states of individual ensemble members that have pronounced three-dimensional (3D) structure, contrasting the homogeneous white noise perturbations applied by Germe et al. (2017).

Ocean initial condition uncertainties and their impacts on climate prediction have been also addressed through bred vectors (Baehr and Piontek 2014) and anomaly transform methods (Romanova and Hense 2015) yielding a weak improvement of prediction reliability at seasonal timescales. Recently, Marini et al. (2016) have achieved a greater ensemble spread for sea surface temperature (SST) during the first 3 years of simulations when oceanic singular vectors are used rather than atmospheric-only perturbations. However, for more integrated measures, such as the North Atlantic SST or the Atlantic Meridional Overturning Circulation (AMOC), the ensemble spread is overestimated initially and decreases over time.

Several studies highlight the strong impact of the 3D structure of ocean state initial errors and emphasize the sensitivity of North Atlantic decadal variability to initial conditions in the deep ocean (Zanna et al. 2011; Palmer and Zanna 2013; Sévellec and Fedorov 2013a, b, 2017). These analyses, based on the singular vectors decomposition (SVD, e.g. Zanna et al. 2011; Palmer and Zanna 2013) or the linear optimal perturbations framework (LOP; Sévellec et al. 2007; Sévellec and Fedorov 2013b, 2017), compute small initial perturbations that induce the maximum response in the system after a specific time. While the SVD method requires solving an eigenvalue problem, the LOP method relies on an optimization problem producing the maximum linear growth of a chosen climatic variable. By construction, both SVD and LOP methods, as applied to the ocean, are based on a linearization of the primitive equations of motion and neglect potential effects of the ocean–atmosphere coupling together with stochastic noise arising from atmospheric synoptic variability. Therefore, assessing the impact of these structures within the full ocean–atmosphere climate system is necessary to better understand their potential for climate prediction.

In this study, we investigate for the first time the impact of LOPs on climate variability in a fully coupled earth system model IPSL-CM5A-LR (Dufresne et al. 2013). We apply the LOP framework maximizing changes in the North Atlantic mean temperature (NAMT) as described in Sévellec and Fedorov (2017). In the ocean model they used, the most efficient LOP induces a NAMT anomaly that reaches its maximum after 10 years. The optimization problem made use of the tangent linear forward and adjoint versions of the ocean component of IPSL-CM5A-LR.

The LOPs dynamics are ultimately related to the excitation of an ocean basin mode identified in the same linear model by Sévellec and Fedorov (2013b). This oscillatory mode involves the westward propagation of subsurface density anomalies across the North Atlantic basin. This propagation impacts the AMOC via thermal wind balance and basin-scale variations of the zonal density gradient. There is evidence of a similar westward propagation in the North Atlantic observations of sea-level height (e.g. Tulloch et al. 2009; Vianna and Menezes 2013), subsurface temperature (Frankcombe et al. 2008), and SST (Feng and Dijkstra 2014) with a comparable basin-crossing time (~ 10 years) as estimated by Sévellec and Fedorov (2013b). It has been also identified in nearly 20 models of the CMIP5 database (Muir and Fedorov 2016). In IPSL-CM5A-LR in particular this oceanic mode exhibits interaction with convective activity, sea ice, and atmospheric circulation (Ortega et al. 2015).

In the present analysis, climate response to the LOP is investigated in terms of changes in NAMT, the AMOC strength, SST, and atmospheric temperature and precipitation. We use ensemble experiments in order to extract the signal of the LOP response from the atmospheric stochastic noise in a perfect model configuration, therefore avoiding pollution of the signal by model drift, and model imperfections. The ensemble experiments, the coupled system and the LOP are described in more detail in Sect. 2. The response of the system to the oceanic perturbations is then described in Sect. 3, while implications for near-term climate prediction are discussed in Sect. 4. Finally concluding remarks are given in the last section.

2 Method

2.1 Model

We use the IPSL-CM5A-LR climate model (Dufresne et al. 2013). It includes the atmospheric general circulation model LMD5A (Hourdin et al. 2013) with a 1.875° × 3.75° horizontal resolution and 39 vertical levels. It is coupled with the oceanic model NEMOv3.2 (Madec 2008) in the ORCA2 configuration corresponding to a nominal resolution of 2°, enhanced over the Arctic and subpolar North Atlantic as well as around the Equator. There are 31 vertical levels for the ocean with the highest resolution in the upper 150 m. It also includes the sea ice model LIM2 (Fichefet and Maqueda 1997) and the biogeochemistry model PISCES (Aumont and Bopp 2006). The coupling between the oceanic and atmospheric components is achieved via OASIS3 (Valcke 2006). The reader is referred to the special issue of Climate Dynamics (volume 40, issue 9–10) for a full discussion of various aspects of this climate model. The characteristics of the oceanic component of the coupled model are also discussed in Mignot et al. (2013).

This model has been used for several decadal prediction studies. In a perfect model context, it exhibits an average predictability limit for the annual AMOC of about 8 years with variations depending on the AMOC initial state (Persechino et al. 2013). The longest potential predictability of SST reaches up to 2 decades and is found in the North Atlantic Ocean. It is related to decadal AMOC fluctuations. These fluctuations are successfully initialized by nudging the SST field to observations (Swingedouw et al. 2013; Ray et al. 2014). This initialization could be further improved, in a perfect model framework, by additionally nudging sea surface salinity (SSS) (Servonnat et al. 2014) and taking into account the mixed layer depth when specifying the amplitude of the restoring coefficients (Ortega et al. 2017). Hindcasts starting from the SST nudged simulations exhibit a prediction skill up to one decade in the extratropical North Atlantic for SST and in the tropical and subtropical North Pacific for the upper-ocean heat content (Mignot et al. 2016).

2.2 General approach

Firstly, we select a 20-year interval (model years 1991–2010) within the 1000-year long pre-industrial control simulation (thereafter CTL) of the IPSL-CM5A-LR model. This specific period is chosen because it does not exhibit strong variability either for the AMOC or NAMT, which both remain within one standard deviation from their 1000-year means. This is necessary to avoid internal variations that may complicate analysing the response to the applied perturbations. Seven ensembles of simulations are conducted using one single starting date—the 1st of January of this time period (model year 1991). All the ensembles are integrated forward for 20 years with a constant pre-industrial external forcing. All ensembles have a random noise disturbance applied to the SST field of the coupler. The applied noise is identical for all ensembles. As this perturbed SST field is only used when SST is passed to the atmosphere during the integration first time step, this perturbation is considered as an atmospheric-only perturbation, as described in Persechino et al. (2013). Germe et al. (2017) showed that this method is equivalent to applying a random white noise to the whole oceanic temperature field. In addition to this atmospheric perturbation, six ensembles utilize full-depth oceanic temperature perturbations. The pattern of these perturbations corresponds to the LOP as computed by Sévellec and Fedorov (2017) using the tangent linear forward and adjoint versions of the same ocean model as in the coupled run. The six ensembles differ only by the magnitude and/or sign of the oceanic perturbation pattern as described below (see Table 1 for details). The seventh ensemble, without any perturbation to the oceanic temperature field, is taken as a benchmark to assess the impact of oceanic perturbations in the other ensembles and will be further referred to as ATM.

Table 1 Summary of the characteristics of the ensembles

Throughout this analysis, the AMOC strength is defined as the maximum value of the annual, zonal-mean streamfunction within 0°N–60°N and 500–2000 m, while NAMT is defined as a full depth average of the annual oceanic temperature over the North Atlantic within 30°N–70°N. The mean state and variability of CTL is assessed from the interannual average and standard deviation for the entire 1000-year time series.

2.3 Oceanic perturbation pattern

The specific pattern of the 3D global oceanic temperature field used to perturb the oceanic initial state of each ensemble was computed by Sévellec and Fedorov (2017) as optimally perturbing NAMT through the LOP methodology. In the analysis they used the adjoint of the tangent linear version of the oceanic component of IPSL-CM5A-LR. More precisely, it was an earlier version of the ocean component (OPA8.2) for which the adjoint version was available at the time of the LOP computation but this difference should not affect the results. This LOP has been rationalized as the efficient stimulation of the least damped oscillatory eigenmode of the tangent linear version of NEMO, fully described in Sévellec and Fedorov (2013a). In particular, its location at depth, away from strong velocities and density gradients (limiting mean- and self-advection, respectively), allows for a long persistence of the anomaly and efficient stimulation of the eigenmode. This eigenmode corresponds to a 24-year oscillatory mode of both the AMOC and the NAMT related to the westward propagation of large-scale temperature anomalies in the North Atlantic. The basin scale propagation influences the AMOC through its impact on the zonal density gradient. Ortega et al. (2015) showed that in the IPSL-CM5A-LR coupled model, this ocean-only mode is maintained by a coupling with a surface mode of variability and potentially excited by the atmosphere. Such coupling allows the intensification of the damped internal mode through the excitation of the deep convection (Sévellec and Fedorov 2015).

By stimulating this internal mode, the LOP provides the most efficient way to generate an anomaly of the NAMT. The LOP pattern depends on the chosen time scale for the transient growth. In this study, we use the LOP maximizing the NAMT response after 14 years in the linear model. In accordance with the lag identified in Sévellec and Fedorov (2013a), corresponding to the time needed for the AMOC to influence the NAMT, we expect an associated maximum response of the AMOC after 8 years. The LOP pattern exhibits the largest magnitudes in the North Atlantic region (Fig. 1), especially in the deep ocean (top versus bottom panels in Fig. 1). These strongest magnitudes of the LOP are furthermore roughly co-located with the areas of strongest temperature variability in the North Atlantic in CTL (black lines in Fig. 1). In Sévellec and Fedorov (2017), both temperature and salinity perturbation patterns are identified. They have a constructive effect on the density anomaly field. In this study, we use only the temperature perturbation as the primary step to understand the response of the coupled system to the LOP. The magnitude of the LOP shown in Fig. 1 leads, in the linear model, to a NAMT response of approximately 43.8 × 10−3 °C after 14 years, which corresponds to roughly one standard deviation of NAMT in the CTL run (not shown). As the LOP is determined by a linear model analysis, its magnitude can be scaled linearly. Scaling factors of 1, 5, 10, 20, − 10 and − 20 are applied, which allow, at the initial date, to sample the whole range of CTL variability in terms of the NAMT index. These six values correspond to the six last ensembles listed in Table 1. The naming of the ensemble reflects this protocol. For example, P20 applies the LOP as shown in Fig. 1, but with its magnitude multiplied by 20 (with a positive sign), while N20 uses a scaling factor of − 20 (negative sign). P01 is therefore the ensemble that uses the LOP exactly as described in Fig. 1, which would lead to one standard deviation response of NAMT after 14 years in the linear ocean-only model (similarly, ATM can be interpreted as LOP0).

Fig. 1
figure 1

The spatial structure of the imposed linear optimal temperature perturbations (LOP, in °C, colour shading) at the ocean surface (top left panel), and at 217 m (top right panel), 1033 m (bottom left panel) and 2768 m (bottom right panel). The amplitudes shown here correspond to the original LOP, i.e. scaled by a factor of 1 (see text for details). Black contours indicate interannual standard deviation of local ocean temperature in the 1000-year long CTL simulation at these depths. The contours are spaced by 0.4 °C within the range from 0.4 to 2 °C at the surface and at 217 m depth, by 0.1 °C from 0.1 to 0.5 °C at 1033 m, and by 0.02 °C from 0.02 °C to 0.12 °C at 2768 m

3 Impact on climate variability

3.1 Ocean response

The climate model ensembles show that the LOP induces a NAMT anomaly reaching its maximum value roughly 10 years later (Fig. 2, top left panels). In qualitative agreement with the adjoint model analysis, it is preceded by a maximum anomaly of the AMOC 5 years earlier (Fig. 2, bottom left and middle panels). The link between these two responses will be detailed below. For both the NAMT and AMOC indices, the magnitude of the response increases linearly with the magnitude of the perturbation (Fig. 2, right panels). The response is significantly different from the ATM ensemble—according to a t test at the 99% confidence level—only for the largest perturbations, i.e. N20 and P20 (Fig. 2, middle and right panels). However, the linearity of the response suggests that the significance for weaker magnitudes can be improved by increasing the ensemble size and therefore the robustness of the statistical test. The AMOC response to the LOP looks slightly asymmetric, being weaker for negative (N10 and N20) than positive (P10 and p20) LOP. However, this asymmetry is not significant when taking into account the confidence interval of the ensemble means (Fig. 2, bottom right panel). Such linearity through the whole range of perturbation magnitudes in a fully coupled ocean–atmosphere system, which includes a significant amount of non-linear processes, is noteworthy.

Fig. 2
figure 2

Left panels: The response of NAMT (top) and AMOC (bottom) to the imposed perturbations for different LOPs’ amplitudes. Colour lines show the time evolution of the ensemble mean for each experiment. The black curve corresponds to the CTL simulation and horizontal dashed lines indicate ± one standard deviation. Middle panels: the same but only for P20 (red line), N20 (blue line) and ATM (grey line), with shading indicating the 99% confidence interval (according to a t test). The vertical black lines in the middle panels indicate the dates that were used to evaluate the magnitude of the response displayed in the right panels. These dates correspond to a 10-year lag for NAMT (top) and a 5-year lag for the AMOC intensity (bottom) after the LOPs were imposed. The time axes refer to model years. Right panels: The magnitude of the NAMT (top) and AMOC (bottom) response after 10 and 5 years, respectively, as a function of the LOPs’ amplitude. Error bars indicate the 99%-level confident interval around the ensemble means. The solid black lines show the best linear fit. Grey shading in the top panel indicates the response magnitude as estimated from the linear ocean model described by Sévellec and Fedorov (2017)

Although linear, the response is also damped by roughly a factor 3 as compared to the response of the linear ocean-only model (Fig. 2, grey shading on the top right panel) and occurs slightly earlier than expected (delay of 10 years instead of 14 years for the NAMT). Such quantitative differences in the response to the LOP in the fully coupled model as compared to the ocean-forced context are expected because of several factors. For example, atmospheric stochastic noise is absent in the oceanic-forced context. Moreover, in the fully coupled model, the perturbation pattern in the surface layer is rapidly distorted and/or damped by air–sea interactions (Germe et al. 2017), which tend to limit the influence of the LOP pattern to its deeper layers. Also, ensemble members differ from each other by their atmospheric states, which may lead to significant differences in air–sea interactions and in the upper ocean. Hence the ensemble average tends to smooth-out the signature of the LOP in the upper ocean. Consistently, the North Atlantic mean temperature of the first 300 m (NAMT300) is very close to the one in ATM during the first 2 and 4 years for P20 and N20 respectively (Fig. 3, top left panel). On the other hand, over the full oceanic depth, NAMT diverges as early as the first year (Fig. 2, top left panel).

Fig. 3
figure 3

Left top panel: Ensemble mean time evolution of NAMT300 for P20 (red), N20 (blue), and ATM (grey) experiments, with shading indicating the 99% confidence interval (according to a t test). The time axis refers to model years. The black curve shows this index for the CTL simulation with black dashed lines indicating ± one standard deviation. The black vertical lines indicate the years selected for the middle and right panels. Left bottom panel: Correlation map between annual T300 at each grid point and the AMOC index in the CTL simulation. Black dots highlight correlations significant at a 95% level. Middle and right panels: Differences in T300 (colours, in °C) between the P20 and ATM ensemble means 1 year (top middle panel), 5 years (top right panel), 10 years (bottom middle panel), and 15 years (bottom right panel) after the LOPs were imposed. The background T300 climatology field in CTL is shown with black contours (contour interval is 2.5 °C). Horizontal red lines highlight the zonal boundaries of the NAMT index (30 and 70°N). Black dots highlight the areas where the plotted ensemble means are different from the ATM ensemble mean at a 95% significance level

Despite the relatively weak initial perturbation in the upper layer, the response of NAMT300 to the LOP is as significant as for the total NAMT (i.e. integrated over the whole water column) after 10 years (Fig. 3, top left panel). Its spatial distribution exhibits a tripole/horseshoe shape (Fig. 3, middle and right panels) that resembles the fingerprint of the response after 5 years to an AMOC intensification in the model (Fig. 3, bottom left panel). This fingerprint pattern is consistent with what can be inferred from SST observations (Dima and Lohmann 2010). This suggests that this upper layer response is mainly driven by the AMOC maximum response to the LOP at 5 years forecast range (Fig. 2 bottom panels). The influence of the LOP on the AMOC has been described by Sévellec and Fedorov (2013b, 2015) in the tangent linear model. As explained above, the involved mode of variability has also been identified by Ortega et al. (2015) in the control simulation using the same climate model (i.e. CTL in this paper). In the present experiments, the LOP imposed in the North Atlantic modulates the meridional density gradient, thereby favouring an acceleration of the AMOC via thermal wind balance. The interaction of the resulting upper-ocean northward flow and the mean meridional temperature gradient gives rise to a temperature anomaly in the upper North Atlantic Ocean. It is the first time that this effect is prognostically tested and highlighted in a fully comprehensive climate model. It confirms the strong sensitivity of the upper ocean to temperature disturbances in the deep ocean, as described in Sévellec and Fedorov (2013a, b, 2017) and validated here in a coupled model. Such impact on the upper ocean suggests some repercussions of the LOP onto the atmosphere in the North Atlantic region.

3.2 Impact on the atmosphere

The impacts of the LOP on the annual mean SST exhibit a tripole pattern (Fig. 4, 1st row) similar to the response of the vertically integrated temperature over the first 300 m (Fig. 3, top middle and right panels). The response to the positive LOP ensemble P20 is stronger and larger scale than its negative equivalent ensemble N20. This is in accordance with the slight asymmetric AMOC response identified in the previous section. The response to the positive LOP ensemble is associated with stronger atmospheric impacts as well (Fig. 4, 2nd to 4th rows). A significant impact is found on the 2-m air temperature (T2M) over the ocean, but also over land in some areas (Fig. 4, 2nd row). Apart from the eastern part of North America, the continental response to the positive and negative LOP is not symmetric. For example, there is a significant response of T2M over the Scandinavia for the P20 ensemble, which is not found significant for N20. A significant impact is found over the western North Africa in N20, while it is found in the eastern North Africa and Middle East regions in P20. These impacts on T2M persist throughout the year but they are stronger in winter than in summer (Fig. 5). In P20, the T2M pattern evolves slightly with the forecasting year, but the warm anomaly in the North Atlantic region persists throughout the first 15 years of the forecasting period (not shown).

Fig. 4
figure 4

Ensemble-mean anomalies for different climate fields in the N20 (left panels) and P20 (right panels) ensembles with respect to ATM, averaged between years 5 to 10 of the forecast. Anomalies are given for annual mean SST in °C (1st row), annual mean T2M in °C (2nd row), summer seasonal mean (June to August) precipitation in kg s−1 m−2 (3rd row) and winter (January to March) sea level pressure in hPa (4th row). Black dots highlight the areas where N20 or P20 ensemble means are different from the ATM ensemble mean at a 95% significance level

Fig. 5
figure 5

Ensemble-mean surface air temperature anomalies (T2M, in °C) for the P20 ensemble with respect to ATM for all years in the 5 to 10-year forecast range. Anomalies are computed for boreal winter (January to March, left panels) and summer (June to August, right panels). Black dots highlight the areas where the P20 and ATM ensemble means are different at a 95% significance level

In accordance with previous finding based on CTL (Persechino et al. 2013), AMOC associated SST anomalies have a significant impact on summer precipitations over the Sahel region (Fig. 4, 3rd row). The positive LOP consistently induces an increase of summer precipitation over the western African Sahel while the negative LOP impacts central and eastern Sahelian region. This asymmetric response is not very surprising considering the asymmetrical SST response. Nevertheless, the details of the teleconnection taking place in the negative case are not fully understood but are beyond the scope of the present study.

Despite these significant impacts on T2M and tropical precipitations, no significant impact could be identified on the major modes of atmospheric variability over the North Atlantic sector, namely the North Atlantic Oscillation (NAO) and the East Atlantic Pattern (not shown). The impact on the winter sea level pressure (SLP) pattern strongly varies with the forecast range and a robust feature of the LOP impacts is difficult to identify at interannual time scales (not shown). When averaging over the 5–10 forecast years, we find a weak, but significant impact over various regions of the North Atlantic (Fig. 4, 4th row). Again, the pattern of the impact differs between the positive and negative LOP. In N20, the pattern has a significant positive anomaly over the Arctic and non-significant negative anomalies over the North Atlantic mid-latitudes, which may be interpreted as a negative NAO-like pattern. The SLP pattern identified for P20 exhibits a zonal dipole opposing the northeastern coast of America with the southeastern European region. This structure does not resemble any well-known large-scale atmospheric circulation pattern from the literature.

4 Discussion: impact on near term climate predictions

In the previous section, it has been shown that the LOP—although computed from the linear version of the oceanic component—successfully excites a subsurface variability mode in the fully coupled system. This mode is known to be associated to subsurface Rossby wave propagation and the associated AMOC enhancement through thermal wind balance. Furthermore, it has been found that the stimulation of this mode has a significant impact on the North Atlantic SST and some atmospheric variables. However, this impact strongly depends on the magnitude of the LOP, going from undetectable signal masked by the atmospheric stochastic noise (e.g. P01, P05) to significant heat anomalies over Europe during several years (P20). In this section, we re-interpret the magnitude of the LOP in relation with the variability of the system, the observational monitoring system in the real world, and a few other ensemble generation strategies in order to give a better insight of the potential usefulness of the LOP for enhancing climate prediction reliability.

4.1 The LOP in the context of IPSL-CM5-LR internal variability

As mentioned in Sect. 2, the magnitudes of the LOP tested in this study sample a large fraction of the NAMT index variability in CTL. This is highlighted in Fig. 6a, where the colour points, indicating the NAMT value for the different magnitudes of the LOP, are over-imposed on the grey shadings that represent respectively one, two, and three standard deviations of NAMT in CTL at interannual timescale. We can see that P01 and P05 magnitudes lie within one standard deviation of the variability from the mean state, which corresponds to very frequent situations, while P20 and N20, on the other hand, rely within two and three standard deviations, and therefore correspond to extreme, and relatively rare events. However, the same analysis, repeated within 4 different oceanic layers (Fig. 6b–e) highlights strong discrepancies within the water column regarding this magnitude. Indeed, all the six LOPs used here, averaged over the first 300 m on the same spatial domain ([30–70°N] in the Atlantic), are very weak compared to the variability of the averaged temperature in the same layer in CTL (Fig. 6b), while they spread over a larger range of the variability in CTL in the deeper layers (Fig. 6d–e). It is at intermediate depth, between 1000 and 2000 m, that the range of LOP magnitudes chosen here is the strongest as compared to the variability of the oceanic temperature in CTL (Fig. 6d). Indeed, within this layer, the LOP strongest magnitude (P20 and N20) is around three standard deviation of CTL. It could therefore be considered as an extreme event: in the assumption of a normal distribution of the NAMT in that specific layer, the probability of such an event would be less than 1%.

Fig. 6
figure 6

a NAMT anomalies and uncertainty ranges, and their contributions from different ocean layers computed for b 0–300 m, c 300–1000 m, d 1000–2000 m, and e below 2000 m. Anomalies are from the LOP experiments (LOP, colour points and black crosses); uncertainty estimates are obtained from the 10-day lagged perturbation patterns (LAG; range in yellow bars), ocean reanalyses (REA; orange bars), and the Argo data (ARGO; light green bars). See text for details. Note that there is no ARGO estimate in e as ARGO floats currently sample the water column only above 2000 m. Grey shadings indicate ± 1, ± 2, and ± 3 interannual standard deviations of the same indices in CTL

Fig. 7
figure 7

The NAMT index based on the Argo float dataset (surface to 2000 m) from a 10-day averaged data and b annual means. Grey shading gives the upper bound on the error based on the study of Taylor (1997). Red shading gives the annual mean estimation of the error when considering time steps within a year as independent

This highlights that the complex 3D pattern of the LOP might create locally very large perturbations as compared to the variability of the system, even though the strongest magnitudes of the LOP are roughly co-located with the strongest temperature variability in the North Atlantic found in CTL (Fig. 1). To investigate the impact of such strong local perturbations, we have generated an additional ensemble, referred to as P20MSK, in which the initial oceanic anomaly is similar to P20 but saturated to 3 standard deviations of the local variability in CTL. The magnitude of the perturbation of this new ensemble in terms of NAMT index is shown in Fig. 6 as a black cross. The perturbation below 2000 m is in particular considerably reduced, although it still reaches 3 standard deviations locally, as in the eastern part of the basin in particular (not shown). In fact, this reduction of the spatial extent of the LOP does not affect significantly the NAMT and AMOC responses (Fig. 2, black crosses in right panels). It indeed stimulates the same Rossby wave propagation mechanism. This suggests that the oceanic response to the LOP is not directly due to its extreme integrated values but rather to its specifically located anomalies.

In summary, the LOPs exhibit a specific 3D pattern, with largest relative magnitudes at intermediate to bottom depths, and a relatively weak perturbation at the surface, when compared to the internal variability. Therefore, while occurrence of such anomalies is very frequent at the surface for all magnitudes that we have tested, their occurrences are extremely rare in the intermediate and deeper ocean. In that respect, P20 and N20 could be seen as extreme events within the North Atlantic Ocean. If a perturbation resembling the LOP was to be detected, one could suspect—although based on this single coupled model analysis—an AMOC anomaly after 5 years, followed by a NAMT anomaly and possible impacts over land. This brings valuable information to assess the North Atlantic climate a few years ahead. It also raises the question of the ability of current monitoring systems to detect such anomalies. This is especially true for the eastern part of the deepest layer (below 2000 m), where the perturbation is very strong, but lies below the maximum depth covered by current Argo floats.

4.2 The LOP in the context of oceanic initial state uncertainties in the real word

Here we compare the LOP to basic estimations of oceanic state uncertainty based on two major data types commonly used to assess the oceanic state and variability: oceanic reanalyses and the Argo float data (Fig. 6, coral and green bars). Our first uncertainty estimation, based on the reanalyses, consists in the integrated (NAMT spatial domain) annual mean temperature differences between GLORYS and ORAS4 (Balmaseda et al. 2013). We chose these reanalyses as they share the same ocean model (i.e. NEMO) as our coupled system therefore facilitating the comparison on similar grids and tools. However, we reckon that this choice likely tends to underestimate the real uncertainties acknowledged in the reanalysis (e.g., Balmaseda et al. 2015; Palmer et al. 2015). The second uncertainty estimation, more directly based on oceanic measurements, uses the 2°-resolution temperature error field of the objective interpolated Argo float dataset described in Desbruyères et al. (2016). Note that to be comparable to the model analysis, both error estimations of the NAMT have been rescaled by CTL variability. The detailed computation of these estimations, and their absolute value (i.e., before rescaling) can be found in appendix 1 and 2. The two estimations give different results, and this already highlights the complexity of assessing oceanic initial state uncertainties and the large uncertainties that remain on these estimations. However, it gives valuable information on the detectability of the LOP.

According to our estimation, in the upper ocean, even for the strongest LOP, magnitudes tested here could not be separated from uncertainty inferred from both reanalyses and Argo data (Fig. 6b). In contrast, at intermediate and deeper layers, the strongest LOPs can be distinguished: below 1000 m, magnitudes of P10 and P20 are larger than the uncertainty of reanalyses and Argo floats (only above 2000 m for the latter). Between 300 and 1000 m, only the largest magnitudes (i.e., P20 and N20) could be distinguished from this uncertainty.

These results have strong implications for climate predictability, the LOP being a source of predictability when it can be detected by the observations. Indeed, in that case, the initial conditions can be correctly assessed in order to phase the subsurface variability mode with the observations, inducing the accurate prediction of its impacts on the surrounding climate. On the other hand, for magnitudes lying below the detectability limit, analysing the LOP’s impact in the climate model may help anticipate uncertainties in climate predictions. These uncertainties could be decreased by extending the monitoring system in the specific regions highlighted by the LOP pattern. In particular, the ocean and the climate were shown to be strongly sensitive to anomalies located below 2000 m, below the current depth of Argo float sampling. This suggests that the deployment of deep Argo floats in the North Atlantic could lead to significant improvements for decadal prediction skills for the North Atlantic region.

Note that the uncertainty estimation done here corresponds to the error on the annual mean oceanic state, while the LOPs correspond to an instantaneous perturbation of the initial state. However, persistence of the LOP can be seen from Fig. 2b: the initial perturbation persists for more than 1 year before generating the anomaly response. Therefore, although it is likely to underestimate the uncertainties on the instantaneous initial state, this comparison still gives useful operational information.

4.3 The LOP for ensemble generation (perturbation) strategies

Taking into account the LOP in the prediction uncertainties can be achieved by perturbing the initial state directly with the LOP to generate an ensemble. However, other perturbation methods might take into account the uncertainty arising from the variability mode associated to the LOP, depending on how the perturbation pattern projects onto the LOP (Sévellec et al. 2017). Random perturbation of the 3D oceanic temperature field arising from white noise local perturbations in each grid box—like used in Germe et al. (2017)—rapidly goes to zero when averaged on a large spatial domain. We have shown in Germe et al. (2017) that this method does not adequately take into account possible deep density structures in the initial state uncertainties. It is thus likely to underestimate the ensemble spread arising from the subsurface variability mode stimulated by the LOP. Another commonly used perturbation strategy of the ocean initial state in near-term climate predictions is based on lagging the oceanic state by a few days (e.g. Hazeleger et al. 2013). We have estimated the magnitude of such perturbations in terms of NAMT using daily time series of the oceanic temperature in CTL. In practice, for each daily oceanic temperature, we compute the difference with the oceanic temperature occurring 10 days before. Then, we compute the NAMT on these instantaneous anomaly fields and take their minimum and maximum values as the range of the initial perturbations arising from this ensemble generation strategy. According to this analysis, the perturbation of the oceanic state due to a 10-day lagged temperature anomaly field is much larger in the surface layer (Fig. 6b, yellow bar) than in the deeper layers where it remains very close to zero, especially below 2000 m (Fig. 6e, yellow bar). This is consistent with the much stronger high frequency variability of the upper ocean. Therefore, the lagging methodology is very unlikely to generate perturbation patterns that project onto the LOP, and thus to excite the subsurface variability mode.

Thus, generating decadal prediction ensemble through LOPs would sample a very different range of initial state uncertainties than other, more traditional methods illustrated in Fig. 6. Practically, this can be achieved by using LOPs of both sign, in addition to atmospheric perturbation for the ensemble generation. In this analysis, the ensemble resulting from merging N10 and P10 exhibits a larger ensemble spread than ATM for the forecast range near the maximal response to the LOP, i.e. 5 and 10 years for the AMOC and NAMT, respectively (not shown). However, this assessment is limited by the fact that the LOP is designed for a specific metric and a specific timescale. Therefore, an ensemble generation based on LOPs as defined in our study is only properly designed to create the largest ensemble spread for the AMOC and NAMT after 5 and 10 years, respectively. This might create an under- or overdispersive predictions regarding other metrics, in particular in other regions, or time scales. This issue is shared with oceanic singular vectors ensemble generation, since the singular vectors also depend on a chosen norm and time scale. In Marini et al. (2016), they found that using oceanic singular vectors gives a better spread for locally assessed metrics during the first year as compared to atmospheric perturbations ensemble generation, while this spread is overestimated for integrated properties such as the AMOC or area-averaged SST. In their analysis, the 3D pattern of singular vectors used to generate the ensemble is not fully described at depth, but their Fig. 3 shows local values of the initial ensemble spread around 0.25 °C in the North Atlantic Ocean at intermediate depth, which is comparable to our local values of interannual standard deviation in CTL. Therefore, prediction uncertainties arising from initial subsurface density uncertainties pattern as identified by the LOP are potentially taken into account by this method.

5 Conclusions

The impact of linear optimal perturbations (LOPs) of the 3D oceanic temperature field for the North Atlantic temperature and for the large-scale meridional overturning circulation has been analysed in a series of model ensembles simulations performed with the IPSL-CM5A-LR climate model. It has been found that the LOPs, as identified in the adjoint version of the tangent linear model of the IPSL-CM5A-LR oceanic component, induce a similar response in terms of oceanic mean temperature and circulation anomalies in the coupled model as in the linear forced ocean model. The response is nevertheless weaker (roughly by a factor 3) and occurs earlier than expected from the linear ocean model. This can be explained by non-linearities in the fully coupled system and the damping terms arising from ocean–atmosphere interactions, which are both absent in the linear ocean model. The computation of LOPs in a fully coupled system would be very challenging. Indeed, it would imply taking into account atmospheric baroclinic and smaller-scale convective instabilities. Within the linear framework used for computing LOPs, such instabilities would not saturate and would dominate the solution, thus contaminating the large-scale ocean response and preventing the determination of the climatically relevant large-scale solutions sought here.

Nevertheless, although the LOPs are based on a linear forced ocean model and have a maximal signature at intermediate depths, they induce a strong SST change and affect atmospheric surface temperature, precipitation and, to a lesser degree, sea level pressure at a 5 to 10-year forecast range. Previous studies based on idealized ocean model configurations have already highlighted the impact of deep oceanic anomalies on the AMOC (e.g. Zanna et al. 2011); however, our study is the first to confirm and quantify such impact in a fully coupled general circulation model. Furthermore, even though our experimental design is rather simple, these results have strong implications in terms of decadal climate predictability. Indeed, they highlight that anomalies in the deep ocean could have significant consequences for the upper ocean and the atmosphere on timescales ranging from interannual to decadal. The fact that the largest amplitudes of the perturbation are found in the deep ocean can be related with the longer persistence of such anomalies in the deeper ocean, where they remain isolated from mean- and self-advection, as well as from damping induced by interactions the ocean mixed layer and the atmosphere. These anomalies persist sufficiently long to maintain the meridional flow and amplify the transient change of the AMOC, which may explain why they are detected as optimal perturbations for this circulation (Sévellec and Fedorov 2015).

The impact of LOPs on ocean heat content is rather linear, whereas the response of SST and atmospheric variables is strongly asymmetric. Regarding the AMOC, its response exhibits a weak asymmetry. Although not significant in our case, this asymmetry has already been observed in the non-linear ocean forced model as a response to SSS optimal perturbations (Sévellec et al. 2008). As explained in Sévellec et al. (2008), this asymmetry may arise from the feedback of density anomalies on ocean vertical mixing. Indeed, a positive (negative) density anomaly will enhance (reduce) vertical mixing and therefore the deep-water formation, resulting in a stronger (weaker) AMOC. Depending on the stratification before perturbation, the positive and negative perturbations may have a different impact that may induce the asymmetry. Here, even though we selected the initial state from a neutral period regarding the NAMT and AMOC variability (cf. section 2), perfect neutrality is nearly impossible to find. Therefore, the asymmetry found in the response might result from the initial state being closer to one sign version of the LOP than the other. Evaluating the impact of particular initial states on the AMOC response would require extensive additional computations. This will be the subject of future work. Likewise, even though an asymmetrical response of the system to LOPs may arise from non-linear feedbacks or more generally from the non-linear interaction of the stimulated linear response with other modes of variability or through non-linear atmospheric and air–sea–ice interaction feedbacks, we cannot reach strong conclusions from our experiments on that aspect.

The SST response to positive LOPs resembles a horseshoe pattern identified in both the IPSL-CM5A-LR model and the observations by Gastineau et al. (2013) as influencing the North Atlantic Oscillation (NAO) in winter. It also resembles the North Atlantic Multidecadal variability (AMV) pattern identified in the same coupled system (Gastineau et al. 2013). The AMV, also known as the Atlantic Multidecadal Oscillation (AMO; Delworth and Mann 2000; Solomon et al. 2011), is known to influence the climate in the North Atlantic region and in particular hurricanes activity (Goldenberg et al. 2001), and precipitations over North America, Europe, and Sahel (Sutton and Hodson 2005; Knight et al. 2006). A large part of its influence over the Euro-Atlantic region seems to be related to its tropical component with a weaker influence of the extratropical SST anomalies (Davini et al. 2015; Peing et al. 2015). However, Gastineau et al. (2016) found a large oceanic influence of the subpolar SST anomaly on the NAO in the IPSL-CM5A-LR model. While the SST pattern associated with the LOP strongly resembles the SST anomaly pattern associated with a negative NAO-like response in Gastineau et al. (2016), we could not identify a clear impact of the LOPs onto the NAO. This could come from a signal to noise ratio issue. Indeed, ensembles of 75-members were used in their analysis, while we are using here 10 members at the most. This highlights the complexity of the influence of the North Atlantic SST on the surrounding climate.

However, our results suggest that density anomalies in the deep North Atlantic could be an oceanic decadal precursor for the AMV and its climatic consequences. This highlights the potential of correct initialization of the full 3D oceanic state to improve climate prediction. Indeed, detecting such anomalies in the real deep ocean could provide a considerable source of predictability for the AMV. Predicting its climatic consequences requires that the AMV impacts, over the continents and in the atmosphere, are correctly represented in the climate models used to perform the predictions. The validity of this latter assumption remains unclear given, for instance the response to an AMV-like pattern is believed to be poorly simulated (Hodson et al. 2009). The upcoming CMIP6/DCPP protocol (Boer et al. 2016) will allow to better evaluate the skill of new generation climate models to represent such teleconnections between the Atlantic SST variations and the atmospheric dynamics. Given the large impacts of the AMV inferred from statistical analysis of the observations, it is possible that a better representation of these teleconnections in future climate models could further enhance the potential climate impact and utility of a precise-enough measurements of deep ocean anomalies.

A comparison of the LOPs with an estimation of the oceanic state uncertainties based on oceanic reanalyses and Argo float data reveals that even the largest magnitudes used here cannot be detected by current monitoring systems in the upper ocean, where the perturbations are the weakest. In contrast, in intermediate and deepest layers, the largest magnitudes (i.e. N20 and P20) stand out of the uncertainty range assessed from observations and reanalyses, suggesting that they could be detected by these systems and therefore initialized in climate predictions.

Finally, our results suggest that a climate prediction starting from an initial state corresponding to an extreme event in terms of the density anomaly in the deep North Atlantic would benefit from using the optimal structure determined in the ocean-only model in the initialization, therefore potentially increasing the prediction skill compared to the average skill in the North Atlantic region. On the other hand, if similar density anomalies were not detected in the observations, they would become a substantial source of uncertainty that needs to be taken into account in climate prediction systems. The best practical way to incorporate these results into decadal prediction experiments will be discussed in future work.