Surface predictor of overturning circulation and heat content change in the subpolar North Atlantic

The Atlantic Meridional Overturning Circulation (AMOC) impacts ocean and atmosphere temperatures on a wide range of temporal and spatial scales. Here we use observational datasets to validate model-based inferences on the usefulness of thermodynamics theory in reconstructing AMOC variability at low frequency, and further build on this reconstruction to provide prediction of the near-future (2019–2022) North Atlantic state. An easily observed surface quantity – the rate of warm to cold transformation of water masses at high latitudes – is found to lead the observed AMOC at 45^∘ N by 5–6 years and to drive its 1993–2010 decline and its ongoing recovery, with suggestive prediction of extreme intensities for the early 2020s. We further demonstrate that AMOC variability drove a bi-decadal warming-to-cooling reversal in the subpolar North Atlantic before triggering a recent return to warming conditions that should prevail at least until 2021. Overall, this mechanistic approach of AMOC variability and its impact on ocean temperature brings new key aspects for understanding and predicting climatic conditions in the North Atlantic and beyond.

1 Introduction

The north-eastward meandering flow of the North Atlantic Current (NAC) dominates the upper-ocean circulation of the northern North Atlantic (Krauss, 1986). It transports relatively warm waters that release heat to the atmosphere as they flow around the Subpolar Gyre (SPG) and the Nordic Seas, ultimately forming North Atlantic Deep Water that propagates in the deep layers via upper and deep western boundary currents (DWBC) and dispersive interior pathways (Bower et al., 2009; Lherminier et al., 2010; see Fig. S1 in the Supplement for domain boundaries and bathymetric features). On top of sequestering physical and biogeochemical properties in the deep seas, the warm-to-cold conversion of water masses and the meridional overturning circulation associated with it drives a significant meridional heat transport. Its variability is thought to be a major cause of temperature and ocean heat content (OHC) shifts in the upper layer of the northern North Atlantic, with important ramifications for ocean–atmosphere interactions and large-scale climate variability (Bryden et al., 2014; Robson et al., 2017). In particular, the most recent reversal of climatic trends in the north Atlantic SPG since 2005 (warming to cooling) has been attributed in numerical models to a decadal weakening of the ocean meridional heat transport across the southern boundary of the SPG (Piecuch et al., 2017; Robson et al., 2016). The recent return of intense ocean-to-atmosphere heat loss (and associated deep convection) since the mid-2010s (Josey et al., 2018; Yashayaev and Loder, 2017) is now suggestive of an ongoing or approaching re-intensification of the circulation and, consequently, a shift to warming condition in the SPG. Overall, the need for a continuous monitoring of the top-to-bottom current field in the SPG has appeared critical to capture the many components of this warm-to-cold transformation. In 2014, international efforts led to the implementation of an in situ mooring array aimed to fulfil such a need – the Overturning in the Subpolar North Atlantic Program (OSNAP; Lozier et al., 2017).

In the commonly used depth space (z), the AMOC_z stream function helps to simplify the complex three-dimensional velocity field of the North Atlantic into a northward flow of about 16 Sv (sverdrup, 1 Sv = 10⁶ m³ s⁻¹) in the upper 0–1000 m or so and a compensating southward flow at depth, connected vertically by the net sinking of surface waters at high latitudes (Buckley and Marshall, 2016; Wunsch, 2002). However, if one is interested in OHC and the dynamics of buoyancy redistribution in the ocean, an estimator of the circulation in density space (σ) must be preferred, which we will note AMOC_σ hereafter. Such an estimator allows us to fully capture transformation of light water masses into denser ones at high latitudes, along both the vertical overturning and horizontal gyre circulations (Lherminier et al., 2010; Pickart and Spall, 2007).

In the absence of diapycnal mixing, the diapycnal volume fluxes associated with the AMOC_σ at a given latitude must relate to air–sea exchanges of buoyancy within isopycnal outcrops north of this latitude (Fig. 1). This thermodynamic balance between the AMOC_σ and its surface-forced component (noted SFOC_σ hereafter), theorized by Walin (1982) and much later verified with numerical models (Grist et al., 2010; Marsh, 2000), suggests key monitoring and predictive skill of AMOC_σ. This was particularly evidenced in low-resolution coupled climate models, which hold a significant lagged relationship between high-latitude surface forcing and overturning circulation at the southern exit of the SPG (Grist et al., 2009). In a follow-up paper, Grist et al. (2014) estimated the surface-forced component of the AMOC in several atmospheric reanalyses and highlighted their overall consistency in the SPG. An independent validation of those surface indices with observation-based time series of the interior circulation is, however, still missing. Moreover, the potential of such proxy-based reconstruction of the AMOC for predicting OHC variability and new climatic reversal in the coming years remains to be shown.

Figure 1 Schematic of the relationship between meridional overturning circulation at latitude ϕ and isopycnal surface σ – the AMOC (ϕ, σ) – and its surface-forced component – the SFOC (ϕ, σ). Arrows show the progressive transformation of waters across increasing density surfaces balanced by buoyancy loss at the air–sea interface and meridional import and export. Inspired from Marsh (2000).

The primary purposes of the present study are (1) to validate with observational data the predictive skill of surface-forced water mass transformation for AMOC variability and (2) to assess the causal link between AMOC variability and decadal OHC changes in the SPG and perform near-future prediction of those quantities. Regional variability will also be documented, with details on the capability of the in situ OSNAP array in monitoring the basin-wide AMOC_σ.

The paper is structured as follows. Section 2 presents the observational datasets and the methodology used to compute AMOC_σ, SFOC_σ, and OHC. Section 3 gathers the main results of the study and Sect. 4 summarizes and discusses them.

2 Materials and methods

2.1 Data

Monthly gridded potential temperature (θ) and practical salinity (S) profiles from four in situ hydrographic datasets were used. Details on those datasets (EN4, CORA, ISHII, and ARMOR3D) are provides in Table S1 in the Supplement. For each product and at each grid point, the θ and S profiles were interpolated to a regular 20 db vertical spacing. Using the TEOS-10 Gibbs-SeaWater (GSW) toolbox, practical salinity was converted to absolute salinity and potential temperature to conservative temperature, and σ₀ and σ₁ (potential density relative to sea-surface and 1000 m, respectively) were computed. Air–sea heat fluxes (radiative and turbulent) and freshwater fluxes (evaporation and precipitation) were obtained from three atmospheric reanalyses (NCEP2, ERA-I, and CERES, see Table S1). Absolute dynamic topography and associated surface meridional geostrophic velocities were obtained from the AVISO platform and combine sea-level anomalies from multi-mission satellite altimeters and mean dynamic topography from GOCE (Gravity field and steady-state Ocean Circulation Explorer), GRACE (Gravity Recovery and Climate Experiment), altimetry, and in situ data (https://www.aviso.altimetry.fr/en/data/products/auxiliary-products/mdt.html, last access: 20 June 2019). Data were downloaded on 12 September 2017 for the years 1993–2015 and in August 2018 for the years 2016–2017.

The various integrated quantities derived from those data products (such as ocean heat content or overturning stream functions – see description below in Sect. 2.2 and 2.3) were then combined into ensemble mean over the period (1993–2017 for altimetry-related quantities, 1985–2017 otherwise), with associated ensemble standard errors computed as $\frac{\mathit{\sigma }}{\sqrt{N-\mathrm{1}}}$, where σ is the standard deviation and N=4 the number of data products used in the mean. This error captures the spread induced by the different methods used as of today to interpolate sparse in situ observations. Further notes on statistical analysis of the reported results (correlation, trend error) are provided in the Supplement.

2.2 Computation of AMOC_σ, MHT_σ, and associated OHC

The 0–2000 m absolute meridional velocities υ at 45^∘ N are derived by referencing in situ estimates of the geostrophic thermal-wind currents with altimetry-derived sea-surface geostrophic velocities, following previously published methodologies (Gourcuff et al., 2011; Mercier et al., 2015; Sarafanov et al., 2012). The latitude 45^∘ N represents the southern geographic boundary of the SPG with the bulk of the light-to-dense transformation associated with the AMOC_σ occurring north of it (see Sect. 3.1). Moreover, the thermohaline fronts (and the resulting relative velocities) at 45^∘ N are relatively well defined in ocean analysis products due to good data coverage, notably at the western boundary. This latitudinal band is therefore chosen as our reference line for computing a realistic AMOC_σ and undertaking its subsequent mechanistic analysis.

The AMOC_σ stream function is obtained by integrating υ zonally and vertically above each σ₁ surface (spaced by δσ=0.025 kg m⁻³). The maximum value of the resulting stream function at the density level σ_m is written as

$$\begin{array}{}\text{(1)}& {\mathrm{AMOC}}_{{\mathit{\sigma }}_m}=\underset{x}{\int }\underset{\mathit{\sigma }<{\mathit{\sigma }}_m}{\int }\mathit{\upsilon }\mathrm{d}x\mathrm{d}z.\end{array}$$

We restrict such a calculation to the 0–2000 m layer as not all products contain data below that depth¹. This threshold is nonetheless deep enough to capture the level of maximum transformation at 45^∘ N, as well as its variability. The AMOC_σ-driven heat transport is estimated as in Mercier et al. (2015):

$$\begin{array}{}\text{(2)}& {\mathrm{MHT}}_{\mathit{\sigma }}={\mathit{\rho }}_{\mathrm{0}}{C}_p{\mathrm{AMOC}}_{{\mathit{\sigma }}_m}\mathrm{\Delta }\mathit{\theta },\end{array}$$

where ρ₀=1025 kg m⁻³, C_p=4000 J kg^−1 ∘C⁻¹, and Δθ is the temperature difference between the upper and lower limbs of the AMOC_σ (i.e. the area-weighted average temperature of water lighter than σ_m minus the area-weighted average temperature of water heavier than σ_m). Note that Δθ was computed from the EN4.2.0 product that provides full-depth temperature profiles. The change in ocean heat content north of 45^∘ N driven by MHT_σ is then estimated as

$$\begin{array}{}\text{(3)}& \mathrm{OHC}(t{)}_{\mathrm{MHT}\mathit{\sigma }}=\underset{t_{\mathrm{0}}}{\overset{t}{\int }}\left(\mathrm{MHT}(t{)}_{\mathit{\sigma }}-\stackrel{\mathrm{‾}}{{\mathrm{MHT}}_{\mathit{\sigma }}}\right)\mathrm{d}t,\end{array}$$

where t is a given year and the overbar refers to a temporal average over the period 1996–2013. This reference period is assumed to represent a climatological equilibrium state around which MHT_σ fluctuates so that positive (negative) anomalies in MHT_σ result in warming (cooling) north of 45^∘ N. As shown in Sect. 3.3, this assumption yields high and significant correlations between OHC_MHTσ and the observed OHC in the SPG.

2.3 Computation of SFOC_σ

The surface-forced component of the overturning stream function SFOC_σ was computed following common practice and methodologies (Marsh, 2000). For each month and each isopycnal σ (spaced by δσ=0.05 kg m⁻³), SFOC_σ was computed as the diapycnal convergence of the diapycnal volume flux driven by surface density flux wherever σ outcrops north of a given coast-to-coast section:

$$\begin{array}{}\text{(4)}& \begin{array}{rl}{\displaystyle }\mathrm{SFOC}\left({\mathit{\sigma }}^{*}\right)=& {\displaystyle \frac{\mathrm{1}}{\mathit{\delta }\mathit{\sigma }}}\int \int \left[-{\displaystyle \frac{\mathit{\alpha }Q}{C_p}}+\mathit{\beta }{\displaystyle \frac{S}{\mathrm{1}-S}}(E-P)\right]\\ {\displaystyle }& {\displaystyle }\mathrm{\Pi }\left(\mathit{\sigma }\right)\mathrm{d}x\mathrm{d}y,\end{array}\end{array}$$

where

$$\begin{array}{}\text{(5)}& \mathrm{\Pi }\left(\mathit{\sigma }\right)=\left\{\begin{array}{cl}\mathrm{1}& \mathrm{for}\mathit{\sigma }-{\displaystyle \frac{\mathit{\delta }\mathit{\sigma }}{\mathrm{2}}}<\mathit{\sigma }<\mathit{\sigma }+{\displaystyle \frac{\mathit{\delta }\mathit{\sigma }}{\mathrm{2}}},\\ \mathrm{0}& \mathrm{elsewhere}.\end{array}\right.\end{array}$$

The quantity within square brackets is the local surface density flux, α and β are the thermal expansion and haline contraction coefficients, respectively, C_p is specific heat capacity of sea water (4000 J kg⁻¹ K⁻¹), Q the net surface heat flux, E the evaporation rate, and P the precipitation rate. Following Marsh (2000), monthly fields of surface temperature (for density computation) and Q are used herein while monthly climatology values for surface salinity S and E−P are used to avoid introducing punctual spurious surface density anomalies due to poor salinity sampling (especially in the early historical record), notably near the continental margins and seasonally ice-covered areas. We note here that the air–sea buoyancy flux in the SPG, and therefore SFOC_σ, is largely controlled by its thermal component (Marsh, 2000). When for a given month (usually during summer), σ does not outcrop north of 45^∘ N, SFOC_σ is set to zero. Annual averages are then obtained for 1985–2017. Even if SFOC_σ is a surface integral statement, maps of transformation rates can be obtained by accumulating the integrand over outcrops (Brambilla et al., 2008; Maze et al., 2009), as shown later in Sect. 3.1.

In order to directly relate SFOC_σ and AMOC_σ, we rely on the assumption of water mass steadiness in the SPG, meaning that the net accumulation of volume within isopycnal layers is considered to be negligible in front of the import of light water to be transformed and the export of dense water after transformation (Marsh, 2000). To verify this hypothesis, we compute $\frac{\mathrm{d}{V}_{\mathit{\sigma }}}{\mathrm{d}t}$, the yearly local change in the volume of discrete isopycnal layer (in Sv), where V_σ is evaluated on 1 January of each year. Averaging this term north of 45^∘ N and summing below the density level of maximum SFOC_σ yields an evaluation of water mass steadiness during each year. As discussed later, this term can be intermittently important but does not dominate the decadal variability so that a direct link emerges between SFOC_σ and AMOC_σ on those relatively long timescales.

3 Results

3.1 The time-mean AMOC_σ and SFOC_σ

The time-mean depth–longitude field of meridional velocity υ at 45^∘ N (Fig. 2a) is dominated by a western boundary current system in good quantitative agreement with direct current estimates carried out in the region (Mertens et al., 2014; Toole et al., 2017). This current system includes the southward-flowing Labrador Current (LC) adjacent to the slope above 800 m (16 Sv), the upper part of the southward-flowing DWBC circa 47^∘ W with increasing velocities with depth (13 Sv), and the surface-intensified northward-flowing NAC (53 Sv) with its recirculation east of 45^∘ W (28 Sv). Meridional velocities are significantly weaker further east in the gyre interior. Integrating zonally the volume transport above discrete σ₁ surfaces yields the (partial) AMOC_σstream function at 45^∘ N, which reaches a time-mean maximum value of 14.3±1.4 Sv at σ₁=32.15 (Fig. 2b; see also Fig. S2A for the AMOC_σ stream function of each individual product). A similar calculation in depth space yields the (partial) AMOC_z stream function at 45^∘ N, which reaches a time-mean maximum value of 9±0.4 Sv at 700 m depth (Fig. S2B). Therefore, about 60 % of the maximum diapycnal volume flux above 2000 m depth at 45^∘ N is associated with a net downwelling in the vertical plane, the remainder being due to dense waters returning at the same depth as that of the inflowing light waters within the horizontal gyre circulation.

Figure 2 Meridional velocity and transport at 45^∘ N (AMOC_σ). (a) Top panel: the 1993–2017 mean longitude–depth velocity field (in m s⁻¹) at 45^∘ N. The σ₁=32.15 isopycnal across which the maximum diapycnal flux occurs is shown in black. Bottom panel: the depth-integrated (0–2000 m) zonally cumulated transport (in Sv) at 45^∘ N with labels as follows. LC (Labrador Current), DWBC (Deep Western Boundary Current), NAC (North Atlantic Current), and RECIRC (NAC recirculation). Shading indicates the ensemble standard error. (b) The mean AMOC_σ stream function at 45^∘ N (in Sv). Shading indicates the ensemble standard error. The dashed blue line at σ₁=32.15 depicts the maximum transformation rate.

The surface-forced component of the AMOC_σ (noted SFOC_σ, see Eq. 4) shows a maximum time-mean value of 15.4±1.8 Sv at σ₀=27.4 (or σ₁≈32), which reflects a light-to-dense flux that primarily occurs along the NAC path in the eastern SPG south of Reykjanes Ridge and to a lesser extent along the western SPG boundary (Labrador Sea), along the Norwegian margins (Fig. 3). This pattern is consistent with recent mooring-based analysis of the diapycnal overturning in the SPG showing a relatively minor contribution of the Labrador Sea to the basin-wide maximum transformation rates (Lozier et al., 2019). This is because the density level of maximum transformation in the Labrador Sea is well below the density level of the basin-wide AMOC_σ (or SFOC_σ). The spatial distribution of the surface-forced diapycnal volume flux within the domain is inferred by evaluating SFOC_σ at two additional key sections: the international Canada–Greenland–Scotland OSNAP and the Greenland–Iceland–Scotland (GIS) sills. The SFOC${}_{\mathit{\sigma }}^{\mathrm{OSNAP}}$ and SFOC${}_{\mathit{\sigma }}^{\mathrm{GIS}}$ stream functions respectively show a maximum transformation rate of 11.2±1.3 Sv at σ₀=27.52 and 5.4±0.4 Sv at σ₀=27.77, in good agreement with independent in situ calculations of the maximum overturning across the OSNAP line and overflow transport estimates at the GIS (Hansen and Østerhus, 2000; Li et al., 2017). Altogether, the three estimates of SFOC_σ across 45^∘ N, OSNAP, and GIS describe the expected decrease in intensity and increase in density of the maximum transformation rate as one progresses northward. We note that the density level of the maximum SFOC_σ at 45^∘ N is slightly lighter than the density level of the maximum AMOC_σ at 45^∘ N. This is because SFOC_σ cannot account for the positive transformation rate due to the entrainment-driven mixing of the subpolar mode waters with the denser overflow waters in the vicinity of the GIS sills. However, the analysis of numerical simulations shows that such a mixing contribution does not largely affect interannual and decadal variability (Marsh et al., 2005), our primary purpose here.

Figure 3 The surface-forced transformation north of 45^∘ N (SFOC_σ). (a) The 1993–2017 time-mean transformation map across the isopycnal surface σ₀=27.4 (in Sv m⁻²), across which the maximum transformation rate north of 45^∘ N occurs. (b) The mean SFOC_σ stream function (in Sv) at 45^∘ N (red) at the OSNAP line (green) and at the GIS sills (yellow). See (a) for section locations. Shading indicates the ensemble standard error. The dashed lines depict the density levels of maximum surface-forced transformation rate north for each domain. As the computation was made using σ₀, the corresponding surface σ₁ values are shown on the right-hand side y axis. The surface integral of the diapycnal volume flux shown in (a) yields the maximum transformation rate through ${\mathit{\sigma }}_{\mathrm{0}}=\mathrm{27.4}:\mathrm{15.4}\pm \mathrm{1.8}$ Sv.

3.2 The variability in AMOC_σ and SFOC_σ

The maximum AMOC_σ time series, displayed as raw and 7-year low-pass-filtered annual anomalies in Fig. 4a (blue lines), shows an apparent 8-year period variability embedded in a linear decrease during 1993–2010 of $-\mathrm{0.24}\pm \mathrm{0.05}$ Sv yr⁻¹ and a subsequent intensification during 2010–2017 of 0.91±0.19 Sv yr⁻¹. Those changes are largely advective (Fig. S4), indicating minor impact of volume (or σ_m) variability on the AMOC intensity. Volume redistribution associated with the formation history of intermediate water masses in the Labrador and Irminger seas can be important but they remain restricted to the lower limb of the AMOC_σ (not shown). We note that the AMOC_z shares a similar variability with AMOC_σ but of weaker amplitude, indicating an important contribution of the horizontal circulation (versus vertical overturning) to the diapycnal volume flux variability at 45^∘ N (Fig. S4). The gyre contribution to AMOC_σ variability at 45^∘ N is also inferred from an independent mooring-based observation of the (400 m bottom) DWBC intensity at 53^∘ N (Zantopp et al., 2017). Although the shortness of the time series (10 years) only allows a suggestive independent validation, the DWBC variability is found to consistently lead the 2004–2010 weakening and the 2010–2014 intensification of the AMOC_σ at 45^∘ N by 3 years (Fig. S4).

Figure 4 The AMOC_σ and SFOC_σ time series. (a) Annual anomalies in the maximum AMOC_σ (blue) and the maximum SFOC_σ (red) at 45^∘ N (in Sv), with the latter shifted 5 years forward (lag of maximum correlation). The reference (time-mean) period is 1996–2013. Thick lines show 7-year low-pass-filtered time series. The right-hand side axis displays the corresponding heat transport anomalies. The original time line for SFOC_σ is given in the top x axis. (b) The 7-year low-pass-filtered time series of anomalies in the maximum SFOC_σ at 45^∘ N (red – shifted 5 years forward), the maximum SFOC_σ at the OSNAP line (green – shifted 4 years forward) decomposed into contributions from the eastern (thin) and western (dashed) basins, and the maximum SFOC_σ at the GIS sills (yellow – shifted 3 years forward). Shading indicates the ensemble standard errors for each variable.

The second independent validation of AMOC_σ bears the mechanistic explanation of its variability. While the maximums of AMOC_σ and SFOC_σ hardly correlate at high frequency, a striking correspondence between their low-pass-filtered variability is found, with the largest correlation obtained when the former lags the latter by 5–6 years (0.94 at the 99 % confidence level), in line with typical advective timescales in the SPG (Bersch et al., 2007) (Fig. 4a – see also the Supplement for details on smoothing and correlation). Therefore, observational data confirm that surface-forced water mass transformation represents a dominant driver as well as an easily derived proxy of low-frequency AMOC_σ changes across the southern exit of the SPG. Departure from an exact match between AMOC_σ and SFOC_σ relates to the influence of the remaining terms in the volume budget equation, namely diapycnal mixing and volume storage within the SPG interior. As shown in Fig. S5, the latter can be non-negligible on interannual timescale but exhibits minor decadal variability.

The 5-year time lag between AMOC_σ and SFOC_σ time series enables prediction of near-future AMOC_σ variability. Here, the low-frequency strengthening of the meridional circulation observed since 2010 is found to continue at a similar rate until 2022, reaching extreme intensities in 2019 and 2020 similar to those observed in the early 1990s. Those extreme events reflect harsh atmospheric winter conditions in the SPG in 2014 and 2015 associated with large ocean-to-atmosphere heat transfer (Josey et al., 2018). As discussed in the next section, this most recent positive trend in AMOC_σ intensity and its predicted persistence until the early 2020s may substantially increase OHC in the SPG in the coming years.

The decadal variability in the maximum SFOC_σ at 45^∘ N has minor contribution from the Nordic Seas and is effectively captured by SFOC${}_{\mathit{\sigma }}^{\mathrm{OSNAP}}$, although the contribution from regions south of the OSNAP line appears important during its most recent intensification since 2010 (Fig. 4b). The variability in SFOC_σ is dominated by changes in the rate of water mass transformation in the eastern SPG basins, in line with recent mooring-based estimates of the AMOC_σ across the OSNAP line (Lozier et al., 2019). The successive 1-year lag between SFOC_σ at 45^∘ N, SFOC_σ at OSNAP, and SFOC_σ at GIS reflects the progressive northward spreading of transformation anomalies across surface of increasing density (see Fig. 3b).

3.3 The variability in OHC and its 5-year prediction

The lagged correlation between surface-forced water mass transformation and the overturning circulation has important ramifications for the monitoring of past, present, and future fluctuations of AMOC_σ, but does not inform on its role in driving decadal OHC variability in the SPG. To infer such a role, the AMOC_σ-driven meridional heat transport at 45^∘ N – noted MHT_σ hereafter – is computed from the maximum AMOC_σ index (Fig. 4a) and the temperature difference between the upper and lower AMOC_σ limbs at 45^∘ N (see Sect. 2.2 and Eq. 2). The time-mean MHT_σ at 45^∘ N during 1993–2017 reaches 0.43±0.04 PW and is balanced by an ocean-to-atmosphere heat transfer of 0.21±0.04 PW, a small long-term change in OHC within the SPG domain of 0.014±0.002 PW, and a northward ocean heat transport across the GIS sills estimated as a residual as 0.20 PW (consistent with independent estimates; Curry et al., 2011; Hansen et al., 2015; Hansen and Østerhus, 2000).

The cumulated anomalies of MHT_σ referenced to the time window 1996–2013 show a high correlation with the observed OHC within the 0–1000 m layer of the SPG (10–70^∘ W, 45–65^∘ N; Fig. 5). In particular, both the 1993–2006 warming and the 2006–2013 cooling of the region are explained well by the contribution of MHT_σ variability at 45^∘ N (r=0.87 at the 99 % confidence level for 1993–2013). This is consistent with previous model-based inferences that the AMOC_σ is a primary driver of decadal temperature changes in the upper SPG (Desbruyères et al., 2015; Grist et al., 2010; Robson et al., 2016). This causal relationship is, however, not verified during 2013–2015, where MHT_σ induces a warming of the SPG, whereas in situ observations indicate that OHC continued to decrease. This apparent discrepancy reflects the strong air–sea heat flux anomaly that drove a sharp cooling of the upper SPG during those years – the so-called “cold blob” (Duchez et al., 2016; Josey et al., 2018). From 2015, atmospheric conditions were back to “normal” and the MHT_σ-driven warming of the SPG could begin.

Figure 5 OHC time series. Detrended anomalies in OHC within the upper SPG (0–1000 m; 10–70^∘ W, 45–65^∘ N; black, in joules) and MHT_σ-driven OHC anomalies north of 45^∘ N (blue, in joules). Shading indicates the ensemble standard errors for each variable. The SFOC_σ-driven OHC prediction for 2017–2022 is shown in red, with its associated error based on the historical predictive skills of SFOC_σ (Eq. 6). The green patch indicates the “cold blob” era driven by extreme air–sea flux events (Josey et al., 2018).

We finally make use of the remarkable 5-year lead of SFOC_σ on AMOC_σ (Fig. 4a) to make a suggestive prediction of AMOC_σ-driven OHC changes between 2017 and 2022. Annually averaged anomalies of SFOC_σ are scaled by the actual interannual variance of AMOC_σ and converted into an anomalous heat transport relative to 1996–2013 with associated OHC anomalies as previously shown. To make the prediction, we simply anchor the resulting 2017–2022 time series to the last observed OHC value of 2017 (red line, Fig. 5). An uncertainty is added to the prediction based on the skill of SFOC_σ in predicting the historical 1993–2017 OHC (red shading in Fig. 5). This uncertainty is the prediction error ∈_lag from the N_lag-year time series, with lag equal to 1 to 5 years:

$$\begin{array}{}\text{(6)}& \begin{array}{rl}{\displaystyle }& {\displaystyle }{\in }_{\mathrm{lag}}=\\ {\displaystyle }& {\displaystyle }\sqrt{{\displaystyle \frac{\mathrm{1}}{N_{\mathrm{lag}}}}\sum {\left\{\left({\mathrm{OHC}}^{y+\mathrm{lag}}-{\mathrm{OHC}}^y\right)-\left({\mathrm{OHC}}_{\mathrm{SFOC}\mathit{\sigma }}^{y+\mathrm{lag}}\right)-\left({\mathrm{OHC}}_{\mathrm{SFOC}\mathit{\sigma }}^y\right)\right\}}^{\mathrm{2}}}.\end{array}\end{array}$$

Owing to the ongoing intensification of the AMOC_σ and its presumed persistence until 2019 or 2020 (Fig. 4a), and under the (hypothesized) absence of extreme air–sea heat flux events in the near future, the present analysis predicts a rapid OHC surge of $\mathrm{1.03}\pm \mathrm{0.57}\times {\mathrm{10}}^{\mathrm{22}}$ J between 2017 and 2021 (Fig. 5).

4 Conclusions

In this paper we have provided observationally based evidence of a tight causal relationship between low-frequency changes in the rate of surface-forced water mass transformation in the eastern SPG, the variability in the overturning circulation at 45^∘ N, and ocean heat content trends in the SPG. The 5-year delay between surface property changes in the SPG and downstream circulation changes suggests good skills for short-term predictability in the region from the sole use of ocean surface and air–sea interface measurements. Here, a strong intensification of the overturning and associated heat transport from 2010 is found to persist until the early 2020s, driving a new significant reversal of climatic condition in the SPG as temperatures rapidly rise from their last minimum of 2017. The extreme winters of 2014 and 2015 appear as key drivers of those recent and upcoming changes in the SPG. They are found to be responsible for rapidly cooling the upper ocean while feeding a 5-year delayed intensification of the overturning through increased light-to-dense transformation, leading eventually to a sharp warming of the domain. We note that the series of oceanic events described herein, from surface-forced water mass transformation to meridional circulation and heat content changes, are only suggestively presented as a forced response to atmospheric variability. Understanding the extent to which they may belong to a more complex loop of coupled ocean–atmosphere interactions is beyond the scope of the present study. Finally, the present analysis confirms the suitability of the international mooring-based OSNAP array for capturing the bulk of interannual and decadal circulation changes driven by air–sea buoyancy exchanges in the whole subpolar area.