North Atlantic Oscillation contributes to the subpolar North Atlantic cooling in the past century

Sea surface temperature (SST) in the subpolar North Atlantic has significantly decreased at a rate of − 0.39 (±0.23\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\pm 0.23$$\end{document}) K/century during 1900–2020, which runs counter to global warming due to anthropogenic forcing. The cooling in the subpolar North Atlantic, known as the North Atlantic cold blob, could be driven by a host of mechanisms involving both the ocean and atmosphere. Here, we present evidence that changes in the atmospheric circulation over the North Atlantic, in particular a centennial trend towards a more positive phase of the North Atlantic Oscillation (NAO), could have contributed to the cold blob. The positive NAO intensifies the surface wind over the subpolar North Atlantic and induces excessive heat loss from the air-sea interface. According to an idealized mixed layer heat balance model, the NAO induced heat loss alone cools the subpolar North Atlantic by 0.26 K/century, which explains 67% of the observed cold blob SST trend. The NAO-induced cooling is partially offset by the warming effect from the East Atlantic Pattern, and the net effect of changes in atmospheric circulation explains 44% of the observed cooling trend. Thus, besides ocean circulation, including the slowdown of the Atlantic Meridional Overturning Circulation, the large-scale atmospheric circulation might have played an equally important role in prompting the century-long SST changes in the subpolar North Atlantic.


Introduction
Since the beginning of the twentieth century, human activities have led to a rising of sea surface temperature (SST) at rate of 0.6 K/century across the globe (Rayner et al. 2003;Junod and Christy 2020).In contrast to this general warming, the SST in the center of the North Atlantic subpolar gyre has been cooling in the past century (− 0.4 K/century as shown in Fig. 1; Hansen et al. 2010;Drijfhout et al. 2012;Li et al. 2022).This outstanding cooling, manifested as a decreasing trend of local SST, is known as the North Atlantic cold blob.The cold blob signifies a unique role of the subpolar North Atlantic in ocean heat uptake (Winton et al. 2013;Marshall et al. 2015;Shi et al. 2018) and has thus generated local and remote climatic impacts.By influencing the baroclinicity of the extratropical atmosphere, the cold blob could modulate North Atlantic storm track activities and the weather patterns downstream in Europe (Woollings et al. 2012;Gervais et al. 2019Gervais et al. , 2020)).In addition, the heat redistribution processes associated with the cold blob have also influenced the location of the Intertropical convergence zone in the Indo-Pacific (Karnauskas et al. 2021) as well as the frequency and intensity of marine heat waves across the Northern Hemisphere (Ren and Liu 2021).Thus, there is a pressing need to address the driving mechanisms of the cold blob.
For decades, the cold blob has usually been hypothesized as an evidence of the Atlantic Meridional Overturning Circulation (AMOC) slowdown in response to increasing greenhouse gases (GHGs) (Rahmstorf et al. 2015;Caesar et al. 2018;Chemke et al. 2020) and decreased aerosol forcing (Ma et al. 2020).The weakening AMOC, in turn, decreases poleward heat transport, and thus cools the subpolar North Atlantic (Delworth et al. 1993).The AMOC slowdown explains the warming deficit in the subpolar North Atlantic projected by climate models under future warming scenarios (Rahmstorf et al. 2015;Menary and Woods 2018;Liu et al. 2020).By then, the AMOC should have weaken significantly due to increased freshwater flux (Sévellec et al. 2017) from sea ice loss (Liu et al. 2019;Liu and Fedorov 2022) and Greenland ice sheet melting (Jungclaus et al. 2006;Hu et al. 2009) and/or changes in ocean temperatures (Marshall et al. 2015;Liu et al. 2017;Gervais et al. 2018;Levang and Schmitt 2020).
However, it has also been argued that the AMOC paradigm may not fully explain the observed cold blob in the past century.One argument is on whether the AMOC has declined during the past century.While temporally short and spatially sparse observations make it difficult to directly assess the long-term evolution of the AMOC, large disagreements and intrinsic uncertainties prevail among indirect estimates from proxy data (Thornalley et al. 2018;Moffa-Sanchez et al. 2019;Rossby et al. 2020;Caesar et al. 2021;Kilbourne et al. 2022).Moreover, the sixth phase of the Coupled Model Intercomparison Project (CMIP6) models simulate a range of historical AMOC trends with the multimodel mean showing a relatively neutral AMOC trend (Weijer et al. 2020).The other argument is that the cold blob might not be a sole result of decreased northward ocean heat transport due to the AMOC decline.A model-based study suggests that the AMOC weakening contributes less than 30% of the cold blob in terms of its intensity and spatial extent (Fan et al. 2021).In addition, the timing of the cold blob precedes the simulated AMOC decline (Drijfhout et al. 2012).The unexplained SST trend by the AMOC has led to alternative explanations of the cold blob.Specifically, changes in wind-driven gyre circulation alter the ocean heat transport in the subpolar North Atlantic and thus contribute to the cold blob (Keil et al. 2020).In addition, changes in the storm activities have induced increasingly more heat loss across the air-sea interface and promoted deep convection simultaneously, which directly and indirectly cool the ocean surface (Li et al. 2022).
These studies indicate that atmospheric forcing might have played a more important role in the subpolar North Atlantic SST than previously thought, especially since atmospheric circulations have undergone substantial changes over the North Atlantic.These changes include an expansion of the Azores High (Li et al. 2012;Cresswell-Clay et al. 2022) and a northward shift of the jet stream and storm track activities (Gulev et al. 2021;Pena-Ortiz et al. 2013;Osman et al. 2021).However, it remains unknown how and to what extent the atmospheric circulation change has contributed to the centennial cooling in the subpolar North Atlantic.In this study, we quantify the direct thermal contribution of atmospheric circulation to the cold blob that is through the heat flux across the air-sea interface.To isolate this atmospheric thermal forcing, we apply an idealized one-dimensional (1-D) mixed layer heat balance model, similar to what Hasselmann (1976) applied to study climate predictability.
The remainder of the text is organized as follows.Section 2 describes the idealized model, analysis methods, and observation-based datasets.Section 3 presents the results regarding changes in atmospheric circulation modes in the North Atlantic, the associated atmospheric forcing field changes, and the SST anomaly (SSTA) trends in the North Atlantic resulting from the atmospheric forcing.Section 4 is discussion, and Sect. 5 is conclusion.

Data, model and methods
In this study, we apply observation-based datasets (SSTA, sea level pressure (SLP), and air-sea heat flux) to quantify the centennial changes in North Atlantic atmospheric circulation modes and the resultant forcing on SSTA.Here, we quantify one atmospheric thermal forcing mechanism that is through the modulation of air-sea heat flux.For this purpose, we first neglect oceanic heat transport terms and simplify the mixed layer heat budget as a 1-D heat balance between the temperature tendency and the net surface heat flux.Thus, the mismatch between our quantification and the observed SSTA trend should be attributed to the processes unresolved in the idealized 1-D model, including but not limited to Ekman transport, gyre circulation, the AMOC, mesoscale and submesoscale eddies, entrainment, and freshwater flux (Buckley and Marshall 2016;Delworth et al. 2017;Su et al. 2018;Zhao et al. 2018;Wills et al. 2019;Oltmanns et al. 2020;Keil et al. 2020;Small et al. 2020;Li et al. 2022).See Sect.2.3 for detailed description of the idealized heat balance model.

Datasets
The five SST datasets used in this study are the 1° × 1° Met Office Hadley Centre's sea ice and SST dataset (HadISST, Rayner et al. 2003), the 1° × 1° Centennial in situ Observation-Based Estimates of SST, version 2 (COBE-SST2, Hirahara et al. 2014), the 2° × 2° Extended Reconstructed Sea Surface Temperature version 4 (ERSSTv4, Huang et al. 2015) and version 5 (ERSSTv5, Huang et al. 2017), and the 5° × 5° Kaplan extended SST (Kaplan et al. 1998;Reynolds and Smith 1994).As in Fig. S1, the spatial pattern and the magnitude of the North Atlantic cold blob differ among the five datasets.For example, ERSSTv4 and v5 show a larger extent of the cold blob than COBE SST does (Fig. S1b, c, d).However, all datasets agree on the cooling trend of annual-mean and JFM-mean SST in the subpolar North Atlantic in contrast to warming elsewhere in the Atlantic basin (Fig. S1).To emphasize the agreement among observational datasets, evenly weighted average of these five datasets (all interpolated to 2° × 2° grids) was used as the best estimation of historical SSTA (Fig. 1).
The monthly mean SLP is from the 1° × 1° ECMWF twentieth century reanalysis (ERA20C, Poli et al. 2016) covering 1900-2010 and the 1° × 1° NOAA/CIRES/DOE 20th Century Reanalysis v3 (20CR, Slivinski et al. 2021) covering 1836-2015.Air-sea heat fluxes are from the 20CR v2 (Compo et al. 2011) which provides heat flux data for 1870-2012, NCEP/NCAR (Kalnay et al. 1996) for 1948-2020, and ERA-5 for 1950-2020.Similar to SSTA datasets, we averaged the three datasets during their overlapping period to generate a best estimation of net surface heat flux ( Q net ) for the past century.Q net is calculated as where Q SW is net down- ward shortwave radiation, Q LW is net upward long wave radiation, and Q SH and Q LE are, respectively, sensible and latent heat fluxes from the ocean to the atmosphere.

Atmospheric circulation modes in the North Atlantic
We apply empirical orthogonal function (EOF) analysis to SLP anomalies during the 1900-2010 period in the North Atlantic sector (20-80°N, 90°W-40°E).The first three EOF modes of annual-mean SLP anomalies are shown in Fig. 2. The leading EOF is the North Atlantic Oscillation (NAO), whose positive phase features a deepened Icelandic Low together with an intensified Azores High (Fig. 2a; Hurrell and Deser 2009).The second dominant mode is the East Atlantic Pattern (EAP), with a center of action south of Greenland (Fig. 2c), which is commonly interpreted as a southward shifted and more zonally oriented NAO pattern (Wallace and Gutzler 1981;Barnston and Livezey 1987).
The third mode is recognized as the Scandinavian (SCA) pattern (Eurasian Type 1 in Barnston and Livezey 1987;Liu et al. 2014), with a center of action near the Nordic Sea (Fig. 2e).These three modes explain over 60% of the atmospheric variability in the cold blob region (Moore et al. 2013).
The modes of variability of January-February-March (JFM) mean SLP anomalies show enhanced centers of action compared to those of the annual mean SLP anomalies, as the extratropical atmospheric circulation is more active dynamically in the boreal winter (Fig. S3).In addition, the two century-long reanalysis datasets analyzed in this study are highly consistent in the EOF modes of North Atlantic SLP in terms of spatial pattern and temporal evolution (Figs.S4, S5 and S6).Even though the centers of the SCA SLP anomalies are located further westward in 20CR compared to ERA20C (Figs. 2e and  S4e) and the principal component (PC) time series correlation is only moderate (0.42; Figs.2f and S4f), the disagreement in   over the Atlantic.The trend is calculated using a linear fit with time that minimizes root mean square error and is based on the averages of five observationbased datasets (HadISST, ERSSTv4, ERSSTv5, COBE-SST2 and Kaplan).The hatched are where more than 80% of the datasets agree on the cooling trend of SSTA SCA modes only marginally affect the quantification of atmospherically forced SSTA trend due to its relatively small mode percentage.As ERA20C and 20CR are consistent with each other in the dominant mode patters of SLP anomaly variation and the key conclusion of this study is independent of the choice of the SLP dataset, and we only present the results based on the ERA20C hereinafter.

Quantify atmospheric forcing on centennial
SSTA trend in the North Atlantic

Ocean mixed layer heat balance
Temperature change in the ocean mixed layer is caused by heat flux across the air-sea interface ( Q net ) as well as vertical and horizontal heat transport mechanisms, including the horizontal heat advection by ocean currents ( Q adv ), entrainment of subsurface water ( Q ent ), and eddy heat diffusion ( Q dif f ).Thus, the mixed layer heat balance can be mathematically formulated as: where T is mixed layer temperature which approximates SST, o = 1024Kgm −3 is seawater density, and is the specific heat of seawater, and H is the monthly climatology of the local mixed layer depth (MLD), which is derived from the EN4.2.1 gridded monthly temperature/salinity profiles (Good et al. 2013).Potential density criteria, i.e., the depth where potential density increases by 0.125Kgm −3 over the surface (5 m) density, are used in the derivation of the MLD (Monterey and Levitus 1997;de Boyer Montegut et al. 2004).It is noteworthy that the year-to-year variation of the MLD is found to have negligible contribution to the simulated SSTA trend based on the idealized model described below (Li et al. 2020(Li et al. , 2022)), primarily due to a trade-off effect between heat inertial of the upper ocean and the effective heat rate atmosphere provides to the mixed layer (Larson et al. 2018;Li et al. 2020).By removing the monthly climatology of each term on Eq. ( 1), the mixed layer temperature anomaly equation is expressed as: To isolate the SSTA change purely due to local atmospheric thermal forcing, we first ignore all the terms involving oceanic heat transport.Resultantly, the temperature anomaly tendency is solely determined by the net surface heat flux anomaly: where T ′ is the mixed layer temperature anomaly, i.e., SSTA.

Decomposition of Q
′ net as a damping and forcing mechanisms on SSTA From the ocean mixed layer heat balance (Eq.3), the net air-sea heat flux anomaly ( Q ′ net ) is a forcing mechanism on SSTA.However, due to its dependence on SSTA, it also serves as a damping mechanism (Stephens et al. 2012).For example, a positive SSTA increases air-sea temperature and humidity differences, which induce a stronger sensible and latent heat flux from the ocean to the atmosphere and thus restore the existing SSTA (i.e., a damping mechanism).The damping and forcing mechanism exerted by Q ′ net can be quantified as: On the right-hand side of Eq. ( 4), the term − T � quantifies the damping mechanism which represents the dependence of Q ′ net on existing SSTA ( T ′ ).The other term Q ′ atmo quantifies the forcing mechanism which is the anomalous heat flux induced purely by atmospheric variables.Both the damping and forcing mechanisms consist of processes involving radiative and turbulent heat fluxes, and thus can be separated as: where turb and rad , respectively, quantify the SSTA damping by turbulent and radiative heat fluxes.
The decomposition of damping mechanisms by turbulent heat fluxes has been described in detail in Li et al. (2020) based on bulk formulae that relate the turbulent heat fluxes to surface wind speed ( ), the air-sea temperature difference ( T − T a ), and the air-sea humidity difference ( q − q a ) as: where a = 1.225Kgm −3 is t he density of air, C D = 1.15 × 10 −3 is the transfer coefficient for sensible and latent heat, (5) and L v = 2.5 × 10 6 JKg −1 is the latent heat of vaporization.
Here, the positive sign of the sensible and latent heat fluxes means from ocean to the atmosphere.Apply the Reynold's decomposition to Eqs. (6, 7) and neglect the higher-order terms, the turbulent heat flux anomalies can be quantified as: Here overbars are the monthly climatology and primes are the deviation from the climatology.In Eq. ( 9), humidity at the ocean surface is saturated and is solely determined by temperature.Thus, q' can be formulated as Eq. ( 8) can be expressed as: Equations ( 8) and ( 10) demonstrate that anomalies in sensible and latent heat fluxes are a function of SSTA (T � ) and atmospheric variables ( . Anomalies in the atmospheric variables may result from internal atmospheric variability or be the response to the underlying SSTA.By assuming a linear relationship between the atmospheric variables and the existing SSTA, the atmospheric response to the SSTA and atmospheric internal variability (residual terms) is separated.As a result, Eqs. ( 8) and ( 10 (10) The SSTA-dependence of Q ′ turb (term in brackets on the right-hand side of Eq. ( 13) provides an SSTA damping mechanism, whose intensity can be quantified by a damping coefficient, turb .As shown above, turb consists of three components: the direct response of turbulent heat fluxes to SSTA ( self ), the response of wind speed to SSTA ( |U| ), and the thermal adjustment of air temperature and humidity to SSTA ( thermal ).The term self is determined by the background wind speed ( can be calculated based on their covariance with SSTA at one month lead time, similar to Frankignoul et al. (1998): With the quantification of turb , the turbulent heat flux anomalies can be partitioned as: . Similarly, the damping and forcing mechanisms by radiative heat flux can be quantified as: . Combining Eqs. 13, 17 and 18, the total damping coefficient is quantified as self + |U| + thermal + rad , that is: According to our quantification, on average is 20Wm −2 K −1 over the North Atlantic and maximizes along the Gulf Stream where high base-state SSTs and strong surface winds promote a more effective SSTA restoration (Fig. 3a).The magnitude and spatial distribution of are largely explained by turb (Fig. 3b), while the damping by radiative processes is an order of magnitude smaller and thus negligible (Fig. 3c).The residual heat flux is independent of the existing SSTA and is determined by the atmosphere, which thus exerts a forcing on SSTA.This atmospheric forcing mechanism via air-sea heat fluxes is quantified as With the decomposition of Q ′ net , Eq. ( 3) is simplified as a stochastically-driven model based on the idea of Hasselmann (1976) applied to assess climate predictability: In our quantification, Q ′ atmo is parameterized as a linear combination of the three EOF modes and a residual function: The simulation is initialized with SSTA in January 1900, and the time resolution is monthly.We parameterize as a normally distributed white noise with its mean equal to zero.The variance of the normal distribution function ( 2) is estimated as the variance of the residual from the linear regression.At each time step, is randomly drawn from the zero-centered normal distribution function ( N 0, 2 ).The ( 17) simulation is repeated 1000 times to quantify the uncertainty range of the simulated SSTA.

Changes to the North Atlantic atmospheric circulation patterns
The first three EOF modes, the NAO (34.1%),EAP (15.2%), and SCA (13.4%), collectively account for 62.7% of the total explained variance of annual-mean SLP anomalies in the past century (Fig. 2).Overall, the spatial patterns of NAO and EAP modes do not show significant seasonal differences, as the annual-mean pattern (Fig. 2) and JFMmean pattern (Fig. S2) show similar centers of action.The SCA mode, however, differs seasonally.The annual-mean SLP anomaly resembles a triple pattern with high pressure (positive SLP anomalies) over Western Europe and North America while low pressure (negative SLP anomalies) over Greenland.The high-pressure center over North America is weakened and even becomes negative in JFM.As a result, the SCA during JFM more resembles a dipole pattern, with the high-pressure center over Western Europe intensifying and expanding toward the eastern subpolar gyre.Accordingly, the low-pressure center shifts westward to the Labrador basin (Fig. S4c).Despite the seasonal differences in SCA mode patterns, the contribution of the SCA to SSTA trend is not significantly affected due mainly to the smaller mode percentage compared to the NAO and EAP.
The atmospheric modes have undergone noticeable changes during the past century.Most importantly, the positive NAO is becoming more dominant, as shown by the significant upward trend of PC1 time series for annual-mean SLP anomalies (1.0 standard deviation (std.) per century, Fig. 2b).Both 20CR v3 and ERA20C datasets agree on a positive trend of the NAO index over 1900-2010, especially a dominance of the positive NAO since the 1940s, consistent with previous studies (e.g., Moore et al. 2013).The tendency from negative NAO towards positive NAO since the 1940s is also shown in the stationbased NAO indices (Hurrell and Deser 2009).In addition, station-based indices and PC-based indices generally agree on the temporal evolution of the NAO in the past century (Fig. S5) and an overall upward trend (Fig. S6).
Moreover, the variance of NAO has increased on decadal to multi-decadal timescales, as shown by running variance maximum occurring in the late twentieth century.Specifically, the maximum 10-year running variance occurs around 1987 (Fig. 4a), and the 30-year running variance reached its maximum (0.06) around the late 1970s (Fig. 4d).The variance increase is even more apparent in the 20CR v3 dataset (Fig. S7).Therefore, the mode of variability of North Atlantic atmospheric circulation is evolving towards a more NAO-like pattern with time.
For the EAP, the PC2 time series displays a neutral trend regardless of seasons, and the annual-mean SLP anomaly shows a 1920-to-1960 dominance of the positive EAP, consistent with other analysis (e.g., Comas-Bru and Hernández 2018; Moore et al. 2013).The variance of the EAP, however, has significantly decreased on decadal timescales, with a maximum before the 1940s (Fig. 4b, d).In contrast, the PC of the SAC exhibits a downward trend (1.6 std.per century; Fig. 2), without changes in the variance (Fig. 4c).
Overall, consistent with previous studies, the past century has witnessed changes in the atmospheric circulation over the North Atlantic, with a more positive NAO, an amplification of NAO mode variance, as well as a reduction of EAP variance (e.g., Hoerling et al. 2004;Hurrell and Deser 2009;Mellado-Cano et al. 2019).The dominance of the positive NAO might have led to more persistent and stronger-thanaverage surface westerlies across the mid-latitude North Atlantic (Thompson et al. 2000;Tamarin-Brodsky and Kaspi 2017), which in turn can promote more heat loss from the ocean surface to the atmosphere, and thus cool the SST.We quantify the thermal contribution of atmospheric circulation changes to the SSTA in the next section.

Atmospheric forcing ( Q atmo ) corresponding to circulation mode of variability
We first quantify the Q atmo field corresponding to the three EOF modes of SLP over the North Atlantic by regressing Q atmo upon the PC time series.Corresponding to the posi- tive NAO phase, the atmospheric field shows a tripole pattern across the North Atlantic, with a weak cooling effect in the tropical North Atlantic, a warming southward of the Gulf Stream, and a significant cooling in the subpolar North Atlantic (Fig. 5a).The Q atmo pattern resembles the tripole SSTA pattern associated with the positive NAO.However, it is noteworthy that the Q atmo is not a result of the NAO- induced SSTA, in that the Q atmo is the residual of Q ′ net with local SSTA's influence removed as damping processes ( − T � ).The Q atmo pattern corresponding to the positive NAO suggests that the observed trend towards a more positive NAO in the past century could lead to a cooling of the subpolar North Atlantic and thus contribute to the formation of the cold blob (Fig. 5a).
The Q atmo corresponding to the EAP and SAC are shown in Fig. 5b and c, respectively.Both SLP modes show the largest loading on Q atmo over the subpolar North Atlantic.The positive phase of the EAP shows a cooling of 6 Wm −2 southward of the Greenland tip (Fig. 5b).In addition, the positive SCA corresponds to a cooling over the Labrador Sea basin and a warming over the Greenland-Norwegian region (Fig. 5c), consistent with a warmer surface temperature during the positive phase of the SCA as identified in previous studies (e.g., Barnston and Livezey 1987;Bueh and Nakamura 2007;Wang and Tan 2020).
For all three modes, the regressed Q atmo pattern is more pronounced in the winter months, as the regression coefficients are about three times the magnitude of the annualmean (Fig. S8a-c).In addition, the tropical signal becomes prominent, suggesting a potential enhancement of tropicalextratropical teleconnection during winter months when the waveguide is relatively southward.
The causes of Q atmo change can be explored by quantify- ing the contribution of surface wind speed ( ), air temperature ( T a ) and relative humidity ( RH ).Mathematically, is calculated using six hourly u, v-wind at surface. is the residual which represents the higher-order processes involving the non-linearity and the covariance between wind and air temperature/humidity in causing Q atmo changes.Further, RH is generally constant near the ocean surface (Dai 2006;O'Gorman and Muller 2010;Schneider et al. 2010), and thus Eq. ( 22) is simplified as: The partial derivative terms are calculated as the linear regression coefficients between the two given variables.
According to the quantification, surface wind speed is primarily responsible for the Q atmo patterns associated with the three atmospheric modes of variability, while the contributions of the near-surface air temperature are regional and an order of magnitude smaller (Fig. 5d-f).Specifically, during a positive NAO, the westerly wind intensifies in the subpolar North Atlantic along with the anomalous easterlies south of 30°N (Fig. 5d), consistent with the established linkage between the NAO and the meridional movement of the North Atlantic jet stream (Woollings et al. 2010;Woollings and Blackburn 2012;(22) Woollings et al. 2018).The resultant wind speed changes induce a tripole pattern in Q atmo , with a 10 Wm −2 cooling in the subpolar North Atlantic, a 4 Wm −2 warming in the subtropical gyre, and a 6 Wm −2 cooling in the tropical oceans (Fig. 5d).The wind-induced Q atmo pattern largely explains the Q atmo during a positive NAO (Fig. 5a), espe- cially over the tropical and subpolar North Atlantic.The minor differences are over the subtropical gyre where the wind-induced Q atmo is about 8 Wm −2 weaker than the total Q atmo change, especially along the extension of the Gulf stream.It is speculated that that higher-order processes, i.e., the covariance between wind and temperature along the storm track (Williams et al. 2007), might also have contributed to Q atmo change in these regions.These higher- order terms are not captured by the decomposition (Eq.23; Fig. 5a, d, and g).
The dominance of wind to Q atmo variability is also mani- fested in the EAP and SCA modes.During a positive phase of EAP, the wind speed changes cool the mid-latitude North Atlantic (~ 45°N) by more than 10 Wm −2 (Fig. 5e), which explains approximately 90% of the total Q atmo change as shown in Fig. 5b.The wind also explains the cooling over the Labrador Sea when a positive phase of the SCA is observed.However, the warming in the Nordic Seas remains unexplained by the wind (Fig. 5c and f), where the provides an approximately 2 Wm −2 warming effect but hardly matches the observed warming (8 Wm −2 ; Fig. 5c).Previous studies have linked the SCA pattern to atmospheric eddy feedback (Bueh and Nakamura 2007), suggesting the covariance between wind and temperature ( ) might be a non-negligible contributor to Q atmo in these cases.Compared to winds, surface air temperature induces smaller Q atmo changes in association with the three domi- nant modes of atmospheric variability (Fig. 5g-i).Overall, the positive phase of the NAO, EAP and SCA is associated with an anomalous cold surface air in the subpolar North Atlantic over their respective center of action (Fig. 5g-i), which provides a cooling effect to SST via Q atmo .Specifi- cally, a positive NAO is associated with 0.6 K cooling in air temperature southward of the Greenland which leads to a −2Wm −2 anomaly in Q atmo (Fig. 5g).The cooling is located over the center of the subpolar gyre with a positive EAP (Fig. 5h), and over the Labrador Sea with a positive SCA (Fig. 5i).Overall, surface air temperature reinforces the cooling over the subpolar North Atlantic induced by winds, even though the temperature-induced cooling is about an order of magnitude smaller than the wind-induced cooling (Fig. 5g-i).
The results based on JFM-mean data (Fig. S8) are consistent with annual mean-based results, except that the magnitude of both wind-and air temperature-induced Q atmo change is approximately three times that of the annual mean counterpart.

North Atlantic SSTA trend in response to Q atmo change
The atmospheric forcing ( Q atmo ) on the North Atlantic SSTA and its contributions to the centennial SSTA trend are quantified using the stochastic model formulated in Eq. 20.
The local forcing field Q atmo is reconstructed using the PCs associated with the three SLP modes (Eq.21; Fig. 5a-c).
According to the model (Eq.20), the combination of the three atmospheric modes has forced a significant tripole pattern of the centennial SSTA trend in the North Atlantic, with a cooling in the southern subtropical gyre, a warming in the northern subtropical gyre, and a most significant cooling in the subpolar gyre and the Labrador Sea (Fig. 5a).The atmosphere forced subpolar SSTA trend differs from observations in terms of the spatial pattern, notably in the Labrador Sea.The Hasselmann model produces a cooling of 0.4 K/century, while the observation demonstrates a warming rate of 1.0 K/ century (Figs.1a and 6a).The mismatch between the observations and the Hasselmann model suggests importance of local and remote ocean dynamics in the Labrador basin, whilst it could also result from the uncertainties in SSTA observations (Fig. 1) and air-sea heat fluxes in this region especially when sea ice is present (Bourassa et al. 2013).Averaged over the Irminger Sea where the observed North Atlantic cold blob is located, the SSTA trend forced by atmospheric circulation modes is − 0.17 [− 0.21, − 0.13] K/century, about 44 [33, 54] % of that observed (Figs.6a and  7).The SSTA trend obtained is consistent with our previous quantification using the total Q atmo trend in the subpolar North Atlantic, in which we showed that changes in Q atmo could have contributed to 52% of the cold blob SSTA trend (Li et al. 2022).Surprisingly, with a highly idealized formulation of mixed layer heat balance (Eqs.3 and 20), our results are qualitatively consistent with the subpolar North Atlantic SSTA trend simulated by a slab ocean model, in which Keil et al. (2020) reported that 30%-40% of the total cold blob SSTA trend is simulated without ocean dynamics.Overall, the simulation by the Hasselmann model suggests that the changes in atmospheric circulation modes in the past century have contributed to the SSTA cooling and the occurrence of the cold blob in the subpolar North Atlantic (Fig. 6a).
Further, we quantify the SSTA trend forced by each individual atmospheric mode as well as their combinations (Fig. 6b-g).With a significant positive trend in the NAO, a cooling in SSTA is induced in the subpolar North Atlantic as expected from the Q atmo (Figs.5a and 6b).With the NAO's dominance as the first EOF mode of SLP anomalies as well as its significant positive trend in the past century (Fig. 2a,  b), the NAO-induced surface cooling explains almost all of the simulated SSTA trend as shown in Fig. 6a.The cooling rate in the Irminger Sea is − 0.26 [− 0.20, − 0.30] K/century, which matches 67 [51, 77] % of the observed cooling trend (Figs. 1a,6b,and 7).This resemblance suggests NAOrelated changes in the air-sea heat flux alone can induce the North Atlantic cold blob during the past century, without explicitly involving oceanic processes.
In contrast, the SSTA trends due to the EAP and SCA are less significant and tend to offset the cooling induced by a more positive NAO.Specifically, the EAP, with an insignificant negative trend of − 0.1 std/century (Fig. 2d), slightly warms the subpolar gyre by 0.07 [0.04, 0.11] K/ century (Figs.6c and 7).At the same time, the more negative SCA (−1.6 std./century) cools the Labrador Sea but warms the eastern subpolar gyre (Fig. 6d).Overall, the combination of the EAP and SCA trend in the past century leads to a 0.09 [0.06, 0.12] K/century warming in the Irminger Sea, counteracting the SSTA cooling forced by the NAO-induced Q atmo (Figs.6g and 7).Such a counteraction limits the south- ern extent of the North Atlantic cold blob.As a comparison, the observed cold blob extends to 40°N, the NAO induced cooling extends to 50°N, while the combination of the three modes leads to a cooling confined to 55°N (Figs. 1, and 6).
Although the observed North Atlantic cold blob extends further south to 40°N, atmospheric forcing principally contributes to the northern portion of the cold blob, as the atmosphere-forced cooling is concentrated on the north of 50°N, which is approximately the latitude of the inter-gyre boundary.SSTAs south of 50°N are predominantly subject to the influence of changes in AMOC-associated ocean heat transport, as suggested by previous model studies (Keil et al. 2020;Fan et al. 2021).Additionally, the southward Ekman transport induced by strong surface westerlies over the intergyre boundary could spread the subpolar cooling anomalies further south.Therefore, we acknowledge the role of winddriven ocean dynamics in the southern extent of the cold blob, which is our ongoing research.While our study suggests a more positive phase of the NAO could contribute to the observed cold blob, our study does not preclude the implicit contribution of the AMOC.The atmosphere and the ocean are tightly coupled in the North Atlantic (Kushnir et al. 2002).Long-term changes in SST are a synergistic effect of natural variability (Delworth et al. 1993;Zhang et al. 2019) and external forcing (Bellomo et al. 2018;Booth et al. 2012), and the impact of atmospheric variability modes is often coupled with the AMOC (Delworth and Zeng 2016;Delworth et al. 2017).For instance, decadal variations in strength of winter convection at Greenland-Iceland-Norwegian Seas are suggested to synchronize with NAO variability (Dickson et al. 1996), which then modulates the North Atlantic thermohaline circulation.A model-based study has found a near-synchronous positive feedback between the AMOC and the NAO, where the atmospheric response to the AMOC is established through AMOC's fingerprint on SST and a surface heat flux damping (Wen et al. 2016).Meanwhile, prolonged positive (negative) NAO could strengthen (weaken) the AMOC by enhancing (inhibiting) deep-water formation through NAO-induced heat flux changes in the subpolar North Atlantic (Delworth and Zeng 2016).In other words, we mainly expect that a more positive NAO would enhance the AMOC and heat transport and counteract the cooling by the direct NAO forcing.On decadal timescales, ocean dynamics is suggested as a crucial bridge between the NAO and subpolar North Atlantic SST variability (Delworth et al. 2017).As such, the atmospheric forcing applied to our mixed layer heat balance model may implicitly include the AMOC-associated component owing to this active AMOC-NAO coupling.

Conclusions
Since the 1900s, the atmospheric circulation in the North Atlantic has experienced a trend towards a more positive NAO with increased variance (Figs. 2 and 4).This positive trend of the NAO index could be due to warming of the tropical Indo-Pacific ocean (Hoerling et al. 2001) and the sea ice loss in the Labrador Sea (Kvamstø et al. 2004;Pedersen et al. 2016).The more dominance of positive NAO leads to stronger cooling effects on the subpolar North Atlantic by promoting heat loss from the ocean surface.The contribution of atmospheric circulation changes to the centennial cooling trend in the North Atlantic cold blob region is quantified based on a stochastic model (Eq.20) that isolates the role of atmospheric thermal forcing via air-sea heat flux ( Q atmo ).According to the quantification, changes in the atmospheric modes of variability in the past century could have contributed a trend of − 0.17 K/century in the Irminger Sea, where the observed cold blob is located (Figs. 6a and 7).In other words, the atmospheric thermal forcing ( Q atmo ) matches 44% of the observed cooling trend, consistent with our previous estimation (Li et al. 2022) and slab ocean model simulations by Keil et al. (2020).Additional sets of controlled simulations suggest that a more positive NAO is a primary contributor to the simulated subpolar North Atlantic cooling (67% of the observed cooling trend), whose cooling effect is marginally offset by the changes in the EAP and SCA.The cooling effects on SST due to the NAO are primarily through the intensification of the jet stream and the wind over the North Atlantic storm track (Fig. 5a).In addition, the reduced surface air temperature reinforces the cooling effects exerted by the stronger surface wind (Fig. 5b).
Our study has pointed out a potential role of atmospheric circulation in forcing the North Atlantic cold blob.However, our results should not be interpreted as a line of evidence against the role of oceanic processes.Atmospheric circulation change induced SSTA trend matches 44 [33,54] % of the observed cooling trend, while the unexplained SSTA trend suggests oceanic processes are indispensable.This result is consistent with Keil et al. (2020) in which the atmosphere drivers only accounts for partial surface cooling in the subploar North Atlantic while a majority of cooling is related to oceanic processes.More importantly, our results do not dispute the role of the AMOC, but rather emphasize the complexity of air-sea interaction in shaping the SSTA trend over the North Atlantic during the past century as the role of atmosphere and ocean differs with respect to time scales and external forcing agents.
It is of caution that the estimate of the long-term AMOC trend could be affected by the relatively short period (decadal or multi-decadal timescales) of AMOC observations/ reconstructions.For example, Fu et al. (2020) combine hydrographic sections, an inverse model and satellite constraints and estimate no overall AMOC weakening since the 1990s.Worthington et al. (2021) apply hydrographic data to an empirical model and estimate no AMOC weakening trend since the year 1981.Jackson et al. (2016) use reanalysis data covering the time period 1989-2015 and suggest the observed AMOC decline at 26.5°N during April 2004-February 2014 as part of AMOC decadal variability.On the other hand, AMOC variability could be complex, as the SST variability in the subpolar North Atlantic likely results from a host of mechanisms (Buckley and Marshall 2016;Foukal and Lozier 2018;O'Reilly et al. 2019;Li et al. 2020), including heat transport by the AMOC (Delworth et al. 1993;Delworth and Mann 2000), stochastic wind forcing (Clement et al. 2015;Cane et al. 2017), air-sea coupling (Brown et al. 2016;Wills et al. 2019), and external forcing (Bellomo et al. 2018;Booth et al. 2012).Sciences (Y.F. and L.L.).R.C. is supported by the National Natural

Fig. 1 a
Fig. 1 a Annual mean and b January-February-March (JFM) mean SSTA trends (shaded; unit: K/century) in the past century (1900-2019) over the Atlantic.The trend is calculated using a linear fit with time that minimizes root mean square error and is based on the averages of five observationbased datasets (HadISST, ERSSTv4, ERSSTv5, COBE-SST2 and Kaplan).The hatched are where more than 80% of the datasets agree on the cooling trend of SSTA

Fig. 2
Fig.2The empirical orthogonal functions (EOFs) of annual mean North Atlantic SLP (upper panels; unit: hPa) and the time series of the corresponding principal components (PCs; lower panel) based on the ERA20C reanalysis.The PCs shown here are scaled by the residuals of the sensible heat flux and latent heat flux, respectively.Adding together Eqs.11 and 12, the turbulent heat flux anomalies are quantified as: and the sensitivity of saturation specific humidity to SSTA, which increases exponentially with background SST according to the Clausius-Clapeyron relation.The terms |U| and thermal depend on the

Fig. 4
Fig. 4 Running variance of PC time series of the a NAO, b EAP, and c SCA as a function of time and running window; d the variance of the three EOF modes with a 30-year running window

Fig. 5
Fig. 5 a-c are Q atmo regressed upon SLP PCs (shaded; unit: Wm −2 ) for the a NAO, b EAP, and c SCA.The Q atmo regression is further decomposed into that cause by (d-f) surface wind speed change( Q atmo  | | | | � ⃗ U | | | |  | | | | � ⃗ U | | | |PC ); and (g-i) surface air temperature (Q atmo T air T air PC ).In a-c, the stippled grid cells are where the regression coefficients are statistically significant at a 0.01 level.The arrows in (d-f) are the anomalous surface wind (vector, unit: ms −1 ) composite on the PCs (PC > 1).The contours in g-i are for the anomalous surface air temperature ( T ′ air ; unit: K) composite on the condition when PC > 1.The bold black contours are zero anomaly isolines and the thin solid (dashed) lines represent T ′ air greater (less) than 0. The contour intervals are 0.2 K ◂

Fig. 6
Fig. 6 SSTA trend (shaded, unit: K/century) forced by the changes in atmospheric modes of variability in the 1900-2017 period based on the simulations by the stochastic model (Eq.20).In the idealized model simulations, atmospheric forcing ( Q atmo ) is parameterized as the combination of the a NAO, EAP and SCA; b NAO; c EAP; d

Fig. 7
Fig.7Box-and-whisker plot of SSTA trend (K/century) averaged over the cold blob region.Each box corresponds to the stochastic model simulations (Eq.20) with Q atmo parameterized as a linear function of atmospheric modes of variability.The boxes represent the interquartile range and the horizontal lines within the boxes represent the median of the simulated SSTA trend.The whiskers end at the maximum and minimum value of all simulations.The gray crosses represent the outliers from the 1000 trails of the randomized simulations.The dashed lines are the 95% confidence interval of the observed cold blob SSTA trend