1
2

Persistent ocean anomalies as a response to northern hemisphere heating3

induced by biomass burning variability

45

Ryohei Yamaguchia, Ji-Eun Kimb,c, Keith B. Rodgersb,c, Karl Steinb,c, Axel Timmermannb,c,6

Sun-Seon Leeb,c, Lei Huangd, Malte F. Stueckere, John T. Fasullof, Gokhan Danabasogluf,

7

Clara Deserf, Jean-Francois Lamarquef, Nan A. Rosenbloomf, Jim Edwardsf

89aResearch Institute for Global Change, Japan Agency for Marine-Earth Science and Technology, Yokosuka,
Japan

10

b

11

c

12
13
14
15
16

Center for Climate Physics, Institute for Basic Science, Busan, South KoreadePusan National University, Busan, South Korea

College of Resource Environment and Tourism, Capital Normal University, China

Department of Oceanography and International Pacific Research Center, School of Ocean and Earth Scienceand Technology, University of Hawai‘i at Mānoa, Honolulu, HI, USAf

Climate and Global Dynamics Laboratory, National Center for Atmospheric Research, Boulder,
CO, USA

17
18
19

Corresponding author: Ryohei Yamaguchi, ryamaguchi@jamstec.go.jp

20

1File generated with AMS Word template 2.0

21ABSTRACT

22

Biomass burning aerosol (BBA) emissions in the Coupled Model Intercomparison Project

23

phase 6 (CMIP6) historical forcing fields have enhanced temporal variability during the years

24

1997–2014, due to the inclusion of satellite records. Recent studies document that the25

corresponding inhomogeneous shortwave forcing can cause changes in permafrost, soil

26

moisture, clouds, which contribute to a net terrestrial northern hemisphere warming, relative to27earlier periods. Here, we investigate the ocean response to the associated hemispherically28

asymmetric warming using two 50-member ensemble simulations conducted with the29

Community Earth System Model version 2 Large Ensemble forced by different BBA emissions

30(CMIP6 default and temporally smoothed over 1990–2020). Differences between the two sub-31

ensemble means show that ocean temperature anomalies occur during periods of high BBA

32variability and subsequently persist over multiple decades. In the North Atlantic, surface

33

warming is efficiently compensated for by decreased northward oceanic heat transport due to34

a slowdown of the Atlantic Meridional Overturning Circulation. In the North Pacific, surface

35

warming is compensated for by an anomalous cross-equatorial cell (CEC) that reduces

36

northward oceanic heat transport. The heat that converges in the South Pacific through the37

anomalous CEC is shunted into the subsurface and contributes to formation of long-lasting

38

ocean temperature anomalies. The anomalous CEC is maintained through latitude-dependent39

contributions from narrow western boundary currents and basin-wide near-surface Ekman

40

transport. These results indicate that interannual variability in forcing fields may significantly

41

change the long time-scale climate background state, presenting a potential uncertainty in

42

CMIP6-class climate projections forced without interannual variability.

43
44

2
File generated with AMS Word template 2.0

45

1. Introduction

46

Aerosol forcings in the latest Coupled Model Intercomparison Project Phase 6 (CMIP6,

47

Eyring et al. 2016) protocol include contributions from biomass burning aerosol (BBA)

48

emissions and a variety of other anthropogenic aerosol emissions (including fluxes from49

industrial activities and land/marine transport, etc.) in addition to natural emissions (dimethyl

50

sulfide, biogenic aerosols, etc.). Although the total amount of global BBA emissions over the51

recent historical period is smaller than other anthropogenic emissions (Fig. 1a), the year-to-

52

year variability of BBA emissions is relatively large. Fires, the main source of BBA emissions,53

have a spatially heterogenous distribution and temporal variations caused by changes in54

vegetation productivity, human activity, climate, and interactions between them (Archibald et55

al. 2013; van der Werf et al. 2017). BBA emissions forcing data in CMIP6 (BB4CMIP6, van56

Marle et al. 2017), spanning 1750 to the present, are based on a multi-model mean flux

57

estimated from the Fire Model Intercomparison Project (FireMIP, Rabin et al. 2017). Satellite-

58

derived observational products derived [Global Fire Emissions Database version 4s (GFED4s),

59

for 1997 onwards], charcoal datasets (in temperate and boreal regions), and visibility records

60(in the Arc of deforestation and Equatorial Asia from the 1960s to 1997) are also incorporated

61

in the BB4CMIP6 (van Marle et al. 2017). The use of observational data which are only

62

available over recent decades, in particular over the satellite era (1997 onwards), provides more

63

globally realistic emissions data but also result in inhomogeneities in the amplitude of temporal

64

variability that are a spurious consequence of patching together disparate types of historical65

reconstructions (Fig. 5 in van Marle et al. 2017). Variability in BBA emissions is greater during

66

the recent period when large-scale observations are more directly incorporated into the67

products (Fig. 1a), although the net amount (long-term mean) of the BBA emissions does not68vary dramatically over time. The spurious change in the amplitude of interannual BBA

69

emissions is a new feature of CMIP6 forcing that was not present in CMIP5, in which decadal

70

mean fluxes were used to construct BBA emissions (Lamarque et al. 2010).

71

3File generated with AMS Word template 2.0

72
73
74
75
76
77
78
79
80
81
82
83

Fig. 1. (a) Globally integrated CMIP6 aerosol emissions forcings from biomass burning(BB, black) and other anthropogenic sources (blue). Temporally smoothed BB emissions overthe period of the large temporal variability (1990–2020, light gray-shaded period) used in thisstudy are superimposed as a red line [BB (SMBB)]. The dark gray-shaded period is when thesatellite (GFED4s) observed data were used (1997–2014). The inserted map shows BB(CMIP6) aerosol emission climatology for 1990–2020 (1011 mol cm2 s-1). (b) Surface(radiative) temperature differences between two sub-ensemble means of a 50-member set
forced by BB (CMIP6) and another 50-member set forced by BB (SMBB) from the Community
Earth System Model 2 Large Ensemble simulations. Statistically insignificant differences at
the 95% confidence level are cross-hatched. See the details of the experimental settings andthe statistical tests in Section 2.

84
4File generated with AMS Word template 2.0

85

Several recent studies have shown that high-amplitude interannual variability in BBA

86emissions can affect not only climate variability for model-simulated fields, but also their87

longer-timescale variability and mean state. This effect is manifested for instance in 2-m air

88

temperature (Fasullo et al. 2022; Kim et al. under review), Arctic sea ice (DeRepentigny et al.

892022), as well as the global hydrological cycle (Heyblom et al. 2022). In particular, in the high-90

latitude Northern Hemisphere (NH), where there are many sources of BBA emissions (Fig. 1a),91

an increase in BBA emissions variability from 1997 onwards causes a net thinning of clouds,

92

thereby enhancing incident shortwave radiation, which in turn warms the surface in the NH

93(Fasullo et al. 2022). The interannual temperature fluctuations generated by the anomalous

94

satellite-era BBA variability can be further rectified into a long-term terrestrial warming

95

through permafrost melting and soil moisture changes (Fig. 1b and Kim et al. under review).

96

As the NH warming expands widely beyond the localized sources of anomalous BBA

97emissions, it induces anomalous Arctic sea ice loss in CMIP6 historical simulations that was98

absent in CMIP5 simulations (DeRepentigny et al. 2022). In turn these processes can further

99

contribute to shifts in the global-scale hydrological cycle (Heyblom et al. 2022).

100

It is generally known that hemispherically asymmetric radiative forcing (provided through101

changes in orbital forcing, albedo, aerosol forcing, or other drivers) leads to interhemispheric102

thermal gradients, and the coupled ocean-atmosphere system works to reduce these gradient103(see the reviews by Chiang and Friedman 2012; Schneider et al. 2014; Kang 2020). In the104

atmosphere, anomalous interhemispheric thermal gradients lead to compensating cross-105

equatorial heat transport (CEHT) by Hadley circulation adjustments and corresponding106

meridional shifts of the Intertropical Convergence Zone (ITCZ) (Broccoli et al. 2006; Kang et107

al. 2008). As for the ocean, compensating CEHT is supplied by anomalies in meridional108

overturning circulation (MOC) and/or anomalies in the shallow subtropical cell structures that109

are partly coupled to perturbations in the atmospheric Hadley Cell (Green and Marshall 2017;

110

Schneider 2017; Kang et al. 2018). A number of studies since the early 1980s have improved111

our understanding of the global response to asymmetric forcing. However the role of the ocean

112

in determining the global response remains relatively unexplored (Kang 2020), despite the fact

113

that its importance is well recognized (e.g., Deser et al. 2015; Yoshimori et al. 2018; Kang et114

al. 2021).

115

The objective of this study is to identify and understand how the anomalous year-to-year

116

BBA CMIP6 forcing translates into anomalous ocean signals that persist beyond interannual
5File generated with AMS Word template 2.0

117

timescales. Our study complements the work of Kim et al. (under review), which focuses on118

the mechanisms that support shorter-timescale rectified land surface responses (warming) over119

the NH. As with the study of Kim et al. (under review), we use output from the Community

120

Earth System Model version 2 Large Ensemble simulations (CESM2-LE, Rodgers et al. 2021).

121

The 100-ensemble-member CESM2-LE is comprised of two sets of 50 ensemble member

122

experiments, integrated from 1850 to 2100 with the CMIP6 Historical (1850–2014) and the123

Shared Socio-economic Pathway 3-7.0 (SSP3-7.0) (2015–2100) scenarios with different BBA

124

forcings. For the BBA forcing, the first 50 members are forced with the default CMIP6 BBA

125

forcing (BB4CMIP6), and the remaining 50 members use temporally smoothed BB4CMIP6

126

for the period of high-amplitude interannual variability (1990–2020) (Fig 1a, details in Section

1272). By taking ensemble mean averages across each 50-member sub-ensemble and differencing

128

the resulting ensemble means, one can deconvolve the forced rectified response and associated129

uncertainty due to background variability. Comparison of the two sets of sub-ensemble

130

averages reveals the presence of ocean anomalies that are induced by the large BBA emissions

131

variability. The ocean temperature anomalies persist over several decades beyond the period

132

of divergent BBA forcing, as well as for the previously reported NH atmospheric warming and133

Arctic sea ice melting (see Section 3 for details).

134

After describing the details of the model experimental design and the statistical methods

135

used in this study in Section 2, we show comparisons of the two sub-ensemble mean fields and136

the atmosphere-ocean responses to the decadal NH warming from the perspective of the energy137

budget in Section 3. In section 4, we provide a mechanistic account of the ocean’s role in138

compensating for the interhemispheric thermal gradient set by the high-amplitude interannual

139

variability of the BBA forcing. The results are summarized along with a discussion of

140

implications in Section 5.

141
142

2. Data and methods

143

We use output from the 100-member LE conducted using CESM2 (Danabasoglu et al.

1442020). CESM2 is a CMIP6-generation Earth System Model developed by the U.S. National

145

Center for Atmospheric Research (NCAR) in collaboration with the broader community. It

146

comprises of the Community Atmosphere Model version 6 (CAM6, Bogenschutz et al. 2018),

147

the Community Land Model version 5 (CLM5, Lawrence et al. 2019), the Parallel Ocean
6
File generated with AMS Word template 2.0

148

Program version 2 (POP2, Smith et al. 2010; Danabasoglu et al. 2012), and the Community Ice

149

Code version 5 (CICE5, Hunke et al. 2015). CAM6 has a resolution of 1.25° in longitude and150

0.9° in latitude with 32 vertical levels. CLM5 has the same horizontal grid as CAM6. POP2

151

and CICE5 share the same horizontal grid, with a uniform spacing of 1.125° in the zonal

152

direction and varying in the meridional direction from 0.65° at high- and mid-latitudes to 0.25°

153

at the Equator, with 60 vertical levels (20 of these are in the upper 200 m of the water column).

154

Initial conditions for the 100-member LE, conducted as part of a collaboration between the155

IBS Center for Climate Physics and the CESM group at NCAR, are a mixture of “macro-” and156“micro-” perturbation initializations (Rodgers et al. 2021). Twenty macro-perturbation initial157

conditions were selected at 10-year intervals from starting points at least 1,000 years into the158

CESM2 pre-industrial simulation (Danabasoglu et al. 2020). By adding small round-off level

159

perturbations to the atmospheric potential temperature field, 20 micro-perturbation initial160

conditions were created from each of the additional four macro-perturbation initial conditions

161(4 x 20 micro-perturbation initial conditions). The four additional macro-initial states were162

chosen based on conditions of the Atlantic Meridional Overturning Circulation (AMOC) (see163

Rodgers et al. 2021 for details). All members were integrated from 1850 to 2100 with common

164

forcings from the CMIP6 Historical forcings for the period from 1850 to 2014 and the SSP3-

165

7.0 scenario for the later period, except for the BBA forcing. Any climate state differences

166

among members originating from differences in the initial conditions are negligible for the167

period of interest in this study (1990–onward, i.e., 140 years after the initialization). Indeed,168

the initial condition dependence of the AMOC, which could have a relatively long memory,169

drops below the detection limit within a few decades (Rodgers et al. 2021, their Fig. S2).

170

The 100-member LE is comprised of two sub-ensembles of 50 members each, with each171

sub-ensemble distinguished by its BBA forcing over 1990–2020 (Fig. 1a). The first 50

172

members are compliant with CMIP6 forcing protocols, thereby including the high-amplitude

173

BBA emissions variability, with this sub-ensemble henceforth referred to as CMIP6. The

174

second group of 50 members uses temporally-smoothed CMIP6 BBA forcing over the period

175

1990–2020 (Fig. 1a) and is henceforth referred to as SMBB. The temporal smoothing was176

performed with an 11-year moving average applied separately for each month of the monthly

177

averaged time series to preserve the seasonal cycle. We emphasize that as a simple moving

178

average was applied, the net integrated BBA emissions over the period of interest are nearly

179

conserved between the two forcings. In addition to the difference in BBA forcing, some slight
7
File generated with AMS Word template 2.0

180

corrections of forcings and code modifications, mainly in CLM5, were applied for the 50

181

SMBB members (see Section 2.4 in Rodgers et al. 2021 for details). The impacts from these182

additional modifications are rather minor with no apparent effects on the climate mean state.

183

For the purposes of this study, we base our interpretation on the assumption that the differences184

between the two sub-ensemble means after year 1990 can be attributed mainly to the185

differences in the temporal variability of BBA between the two sub-ensembles over 1990–2020.

186

To determine if a difference between sub-ensemble means is statistically significant, we187

use a bootstrap-like method to utilize information on background internal variability which can

188

be estimated from the LE members. We first randomly select two groups of 50 members from189

the 100-member LE and compute the difference between the means of the two groups.

190

Repeating the computation 10,000 times, we calculate the standard deviation of the 10,000

191

differences, with the mean converging to zero as the number of trials increases. For the case

192

where the difference between the CMIP6 50-member mean and SMBB 50-member mean

193

(referred to as the 50-50 difference) is outside the range of twice the standard deviation, the194

50-50 difference is defined as statistically significant at the 95% confidence level, i.e., it is not

195

likely to occur by chance due to the background internal variability.

196
197

3. Impacts of high-amplitude interannual variability in BBA emissions

198

a. Anomalous heat entered into the Earth System

199

8File generated with AMS Word template 2.0

200
201
202
203
204
205
206
207
208
209
210
211

Fig. 2. The difference in sub-ensemble means between the 50-member CMIP6 and 50member SMBB simulations for (a) globally integrated top of atmosphere radiation imbalance(TOA, blue line) and ocean heat uptake (OHU i.e., net ocean surface heat flux, black line), (b)interhemispheric difference in column-mean atmospheric temperature (Northern Hemisphereminus Southern Hemisphere), (c) Northern Hemisphere mean atmospheric temperature profiles,(d) Arctic sea ice (SI) volume, (e) global mean steric sea level change (SLC, estimated fromsea water density changes), (f) AMOC maximum transport at 40°N, and (g) global mean ocean
temperature profiles. In (a), the cumulative sum of TOA difference from 1989 (red line withright axis) is also shown. The time interval surrounded by two green lines corresponds to thetime interval where CMIP6 and SMBB forcings differ. The estimated range of standard
deviations (σ) based on the bootstrapping analysis are shown in (a), (b), and (d–f) (light gray
9File generated with AMS Word template 2.0

212
213

shadings for ±σ and dark gray shadings for ±2σ). Statistically insignificant differences at
the 95% confidence level are cross-hatched in (c) and (g).

214
215

Statistically significant CMIP6 vs. SMBB differences in several climatological metrics are216

seen during and/or beyond the period of different BBA forcing (1990–2020) in Fig. 2. In terms

217of the global energy budget, a relative excess of heat enters the Earth system in the form of218

radiative fluxes through the top of the atmosphere (TOA, up to 0.4 W m-2 as an annual mean)219

in the CMIP6 members (blue line in Fig. 2a). The total amount of heat accumulated from the220

beginning of the period with different forcing reaches up to 20 ZJ (red line in Fig. 2a). The

221

anomalous heat input from the top of the atmosphere enters as an anomalous shortwave

222

radiative flux, mainly within the NH high latitudes (Figs. 3a–3c, and Fig. 4a). The heat that223

enters the Earth system serves to immediately warm the atmosphere (Fasullo et al. 2022) and224

to melt Arctic sea ice (DeRepentigny et al. 2022) through the rectified decadal NH warming

225(Kim et al. under review) (Figs. 2c and 2d). In contrast to the rapidity with which the heat enters

226

the system, it is radiated back to space rather slowly (red line in Fig. 2a). Whereas the heat227

largely enters the Earth system in the NH high latitudes, with a maximum spanning from228

Siberia to North America (Fig. 4a), the regions of anomalous emissions back to space as

229

longwave radiation over 2015–2030 occur principally on the eastern flanks of the North Pacific

230

and North Atlantic Oceans (Fig. 4b).

231

10File generated with AMS Word template 2.0

232
233
234
235
236
237
238
239
240
241
242

Fig. 3. The difference in CMIP6 (50 members) and SMBB (50 members) ensemble meansof zonally-integrated (a–c) top of atmosphere radiation imbalance (TOA) and (d–i) ocean heatuptake (OHU). TOA consists of shortwave radiation (b) and longwave radiation (c)components and the modeled OHU is decomposed into components of shortwave radiation (e),longwave radiation (f), latent heat of evaporation (g), sensible heat (h), and other heat from sea
ice process and river runoff (i). 11-year running mean is applied for all panels. Positive
differences indicate anomalous downward flux (TOA: to the atmosphere, OHU: to the ocean)in CMIP6 members. The time interval between the two green lines corresponds to when CMIP6
and SMBB forcings are different (1990–2020). Statistically insignificant signals at the 95%
confidence level are cross-hatched.

243

11
File generated with AMS Word template 2.0

244
245
246
247
248
249
250
251
252

Fig. 4. Differences in top of atmosphere radiation imbalance (TOA) and ocean heat uptake(OHU) fluxes between the CMIP6 and SMBB sub-ensemble means. Spatial patterns are shownfor the periods of absorbing heat (2000–2010, a and c) and releasing heat (2015–2030, b andd), defined based on the global mean TOA flux differences (Fig. 2a). Positive differencesindicate anomalous downward fluxes (TOA: to the atmosphere, OHU: to the ocean) in CMIP6
relative SMBB. 100-member ensemble mean climatologies of sea surface height [SSH (cm)]for each period are superimposed in (c) and (d). Statistically insignificant differences at the
95% confidence level are cross-hatched.

253
254

Less than 10% of the total anomalous heat input during 1990–2020 goes to warm the255

atmosphere and melt sea ice. The vast majority of the anomalous heat (~92%) is absorbed by256

the oceans (black line in Fig 2a). Preferential ocean heat uptake (OHU, i.e., net air-sea surface

257

heat flux) regions include the band represented by the NH trade winds (near 20°N), winter deep258

convection regions near the western boundary currents (WBCs) in both the North Pacific and259

Atlantic, and along the Arctic sea ice margin (around 60°N) where sea ice retreats and hence260

heat exchange is enhanced (Figs. 3d and 4c). Shortwave radiative flux anomalies are the main

261

contributors to anomalous ocean heat uptake (Figs. 3d–3i) in the NH high-latitude. The heat262

taken up by the ocean raises the global mean ocean temperature between the surface and ~700263

m depth (Fig. 2g) and expands the ocean, as can be seen with the global mean steric sea level
12File generated with AMS Word template 2.0

264change of 3 mm (Fig. 2e). Previously reported atmospheric and cryosphere anomalies induced

265

by the high-amplitude BBA emissions (Fasullo et al. 2022; DeRepentigny et al. 2022; Kim et266

al. under review) disappear soon after the end of the period where CMIP6 and SMBB forcings267

are divergent (e.g., Figs. 2c and 2d). In contrast, oceanic anomalies remain in the system for at

268

least several decades, acting as a longer-term memory in the system (e.g., Figs. 2e–2g). In the269

next section, we investigate the mechanisms responsible for generating these long-lasting270

ocean anomalies by tracking ocean heat anomalies within the ocean interior.

271
272

b. Change in atmospheric and oceanic interhemispheric heat transport

273

Increases in the interhemispheric thermal gradient (Fig. 2b) due to hemispherically274

asymmetric warming anomalies lead to coordinated adjustments of the atmosphere and ocean

275

meridional heat transports (MHTs) that act to reduce the gradient, as has been documented in276

previous studies (e.g., Deser et al. 2015; Green and Marshall 2017; Kang et al. 2021). For our277

case, where the ocean transports heat northward and the atmosphere transports heat southward278

across the equator under climatological conditions (Fig. 5a, inlay figure), the NH decadal279

warming anomaly induced by high-amplitude variability in BBA emissions reduces the280

northward ocean heat transport and increases the southward atmospheric heat transport at the281

equator (Fig. 5b). As a result, the excess heat added to the NH is compensated or transported282

to the Southern Hemisphere (SH), thus contributing to decreasing the interhemispheric thermal

283

gradient. Changes in the atmospheric CEHT are linked to changes in the Hadley circulation

284(Kang et al. 2008; Schneider et al. 2014), which is further manifested by the meridional shift285

of the ITCZ. Also in the present case, a northward shift of the ITCZ is clearly identified in the286

Atlantic and Indian Ocean sectors (Fig. S1b). Changes in the wind-driven component of the287

ocean CEHT are coupled via the trade winds to changes in the Hadley circulation (Green and288

Marshall 2017), and for this case the ocean CEHT anomaly occurs in the same direction

289(southward) as those in the atmosphere.

290

13
File generated with AMS Word template 2.0

291
292
293
294
295
296
297
298

Fig. 5. (a) 100-member ensemble mean climatology (2000–2020) of the meridional heattransport (MHT, positive northward) by the atmosphere (ATM) and ocean (OCN). The ocean
MHT is decomposed into the Eulerian mean circulation component (OCN_EM) and othercomponents (OCN_Other) including the mesoscale eddy (bolus) transport, horizontal diffusion,and parameterized sub-mesoscale transport. The Eulerian component in the Atlantic(OCN_EM_At) is shown separately from the Eulerian component over the global domain(OCN_EM_Gl). The inserted panel in (a) shows an enlarged view of the climatological MHTs
14
File generated with AMS Word template 2.0

299
300
301
302
303
304

near the equator. (b) Differences in cross-equatorial heat transport (CEHT) between CMIP6
and SMBB sub-ensemble means. Note that the OCN and OCN_EM_Gl lines in (b) almost
overlap each other. An 11-year running mean is applied to the time series. The gray shaded
period corresponds to when CMIP6 and SMBB forcings are different (1990–2020). The
estimated two standard deviations ranges (±2σ) based on the bootstrapping are shown asdashed lines.

305
306

Responses of the oceanic and atmospheric CEHT to interhemispheric forcings are sensitive307

to the geographical location and intensity of the forcings (Hawcroft et al. 2018; White et al.

308

2018; Yu and Pritchard 2019). In the present case, the magnitude of the CEHT anomaly is309

approximately five times larger in the ocean than that in the atmosphere, with no phase

310

difference between the two (Fig. 5b). The southward ocean CEHT anomalies are mostly

311

induced by the contribution from the Eulerian mean (EM) currents, and the other components

312[mesoscale eddy (bolus) transport, horizontal diffusion, and parameterized sub-mesoscale313

transport] do not contribute strongly. For the climatological mean behavior of CESM2, CEHT

314

for the Atlantic Ocean is positive (northward heat transport) and CEHT for the Indo-Pacific315

domain (OC_EM_Gl minus OCN_EM_At in Fig. 5a) is negative (southward heat transport).

316

However, the directions of the BBA-induced CEHT anomalies are southward in both ocean

317

basins. Of the southward ocean CEHT anomaly, 40% occurs in the Atlantic, and the remaining

318

60% in the Indo-Pacific (Fig. 5b), with this result standing in contrast to a previous study that319

emphasized the Atlantic contribution under the NH cooling (Yu and Pritchard 2019).

320
321

4. Mechanisms for the long-lasting ocean anomalies

322

a. Fate of heat taken up by ocean

323

The fate of the heat absorbed by the ocean is identified by applying a budget analysis for

324

ocean heat content (OHC) (Fig. 6). For each basin, changes in ocean heat content (dOHC/dt,

325

green lines in Fig. 6) should be balanced by the sum of the air-sea heat flux at surface (SHF,

326

orange lines in Fig. 6) and the ocean heat transport convergence (OHTC, blue lines in Fig. 6).

327

It should be noted that the OHTCs are estimated directly from temperature and velocity fields

328

only for those basins where the relevant model output was saved (only for the North Pacific,

329

North Atlantic, Southern Ocean, and Low-latitude Atlantic); otherwise, they are estimated as a

330

residual in the budget analysis. For the full global ocean domain (Fig. 6a, no OHTC term exists

331

by definition), during the period of the high-amplitude variability in BBA emissions (1990–
15
File generated with AMS Word template 2.03322020), ocean heating of ~6 x10-2 PW occurs, resulting in up to 18 ZJ of anomalous OHC in the333

CMIP6 members. Subsequently, the global OHC anomaly decreases more slowly than its

334

increase. Although the pattern of ocean heat loss is less coherent than it is for ocean heat uptake,

335

relatively strong heat release can be identified along the Kuroshio/Oyashio extension region336

and in the eastern Tropical and South Pacific (Fig. 4d). These are due to decreases in incoming

337

shortwave radiation and anomalous latent heat fluxes, respectively (Figs. 3d–3i). While the338

ocean steadily loses heat after the end of the period, the global OHC anomaly remains

339

statistically significant at the 95% confidence level until the end of the simulation (year 2100)

340

in association with accumulation of heat below 300 m depth (Fig. 6a).

341

342
343
344
345
346
347
348
349

Fig. 6. Budgets for the ocean heat content differences between the sub-ensemble means of
CMIP6 and SMBB members. For (a) the global ocean and for (b–h) each ocean basin (its spatial
map is inserted as the upper central panel), changes in the 50-50 differences in the full-depth
ocean heat content (dOHC/dt), the sea surface heat flux (SHF), and the ocean heat transport
convergence (OHTC) are shown as green, orange, and blue lines with the left y-axis,
respectively. The residuals (Res.) are also shown as gray lines with the left y-axis if available.
Black solid lines are 50-50 differences in the full-depth OHC (with the right y-axis). In the16
File generated with AMS Word template 2.0

350
351
352
353
354
355
356
357

panels for global (a), Pacific (d and g), and Atlantic (e and h), 50-50 differences in OHCs above
300 m and below 300 m are also shown as dotted and dashed lines with the right y-axis,
respectively. The Southern Ocean is defined as oceans south of 34°S, and the Pacific andAtlantic are separated into their respective northern and low-latitude parts at 20°N. An 11-year
running mean is applied to all timeseries. The white dots on the black, blue, and orange lines
indicate statistically significant differences at the 95% confidence level based on thebootstrapping method. The gray-shaded period corresponds to the time interval where CMIP6
and SMBB forcing are different (1990–2020).

358
359

The North Pacific (Fig. 6d) and North Atlantic (Fig. 6e), where anomalous surface heating360

due to anomalous BBA forcing is the largest, show seemingly similar responses. Surface

361

heating in the North Pacific and North Atlantic reaches up to 2x10-2 PW during the time of362

high-amplitude BBA emissions variability, resulting in a corresponding local increase in OHC

363(~ approximately year 2016). In both regions, the negative OHTC anomalies (divergence of364

ocean heat transport) following the anomalous positive SHF perturbations dominate the365

dOHC/dt. The OHC anomalies themselves return to a statistically insignificant level by the366

year 2020 for the North Atlantic and by the year 2030 for the North Pacific, although the367

negative SHF anomalies also partially contribute to the OHC decrease in the North Pacific.

368

Contemporaneously with the negative OHTC anomalies in the North Atlantic (~ year 2040),

369

the OHC in the Arctic Ocean is significantly increased by the positive OHTC anomaly (Fig.

370

6b). On the other hand, positive OHTC anomalies in the low-latitude Atlantic (~ year 2020)

371

are subsequently offset by the negative SHF anomalies, and therefore OHC anomalies do not372

increase to a statistically significant level (Fig. 6h). Positive OHTC anomalies in the low-373

latitude Pacific, which occur at the same time as negative OHTC anomalies in the North Pacific374(in approximately years 2000–2010), cause positive OHC anomalies which are the only

375

statistically significant OHC anomalies among ocean regions at the end of the simulation (Fig.

376

6g). The long-lasting low-latitude Pacific OHC anomalies are largely confined to the377

subsurface (deeper than 300 m) and account for the largest fraction (50%) of the remaining378

global OHC anomalies at the end of the simulation.

379

b. Ocean basin response for compensating and/or sequestering the heat

380

Although for both the North Pacific and the North Atlantic the anomalous heating381

associated with higher-amplitude variability BBA emissions is compensated by the ensuing

382

OHTC anomalies, the meridional extent of the ocean heat transport anomalies is notably

383

different between the two ocean basins. During the latter half of the high-amplitude BBA

384

emission period (2000–2020) when the OHTC anomalies are large, the global ocean meridional
17
File generated with AMS Word template 2.0385

heat transport (MHT) anomalies extend broadly from 60°N to the south of 30°S with the386dominant contribution from the Eulerian mean component to its southward anomalies (Fig. 7a).

387

The MHT anomalies for the North Atlantic Ocean are responsible for a portion of the global

388

ocean MHT anomalies broadly from 60°N to at least 30ºS. On the other hand, the negative

389

MHT anomalies in the Indo-Pacific basin, mainly caused by the Eulerian mean component (the390difference between OCN_EM_Gl and OCN_EM_At in Fig. 7a), are confined within 30°S–

391

30°N.

392

The dominant Eulerian mean components in the MHT anomalies in each basin correspond393

well to changes in the MOC over each basin. The meridionally-uniform Atlantic MHT

394

anomalies share their spatial characteristics with the weakening of the AMOC (Fig. 7b). On395

the other hand, in the Pacific the shallow Subtropical Cells (STCs, McCreary and Lu 1994) at

396

low latitudes (30°S–30°N) weaken in the NH and strengthen in the SH during the latter half of397

the period of the high-amplitude BBA emissions (Fig. 7c). The combined effect of the STC

398

anomalies in both hemispheres can form an anomalous Cross-Equatorial Cell (CEC) which is399

comprised by ascending (descending) branches in the NH (SH) and southward (northward)

400

currents in surface (subsurface) waters. Southward surface current anomalies at the equator,

401

corresponding to the upper branch of the CEC, are consistent with the southward MHT

402

anomalies in the Pacific Ocean. These basin contrasts of the MOC responses has been identified403

as well in previous studies that evaluated the response of climate models to hemispherically404

asymmetric cooling (e.g., Yu and Pritchard 2019; Kang et al. 2021)

405

18
File generated with AMS Word template 2.0406
19File generated with AMS Word template 2.0

407
408
409
410
411
412

Fig. 7. Differences in sub-ensemble means between CMIP6 and SMBB for (a) meridional
heat transport (MHT) and (b and c) meridional overturning circulation (MOC, positive
clockwise). Line colors in (a) are the same as in Fig. 5. 100-member ensemble mean
climatology of MOC streamlines (2000–2020, Sv) is superimposed as contours in (b) and (c).
Statistically insignificant differences at the 95% confidence level are cross-hatched in (b) and
(c).

413
414

1) ATLANTIC

415

From a temporal perspective as well, the North Atlantic OHTC anomalies in Fig 6e are416

consistent with the AMOC anomalies. In the sub-ensemble members forced by the high-417

amplitude variability in BBA emissions (CMIP6 members), the AMOC transport, assessed at418

the latitude where it is maximum (40°N), is ~0.7 Sv smaller than that in the SMBB sub-419

ensemble both during the period where CMIP6 and SMBB forcing differ and thereafter (Fig.

4202f). The AMOC anomalies are also largest at 40°N, and extend meridionally across the Atlantic

421

basin from the 2000s to the 2030s, and eventually disappear completely by the 2040s (Fig. S2).

422

In terms of the heat budget, the AMOC anomalies compensate for the anomalous surface heat423

input into the North Atlantic by reducing mean northward ocean heat transport (i.e., a424southward ocean MHT anomaly is induced).

425

The weakening of the AMOC during and after the period of the high-amplitude variability

426

BBA emissions is likely due to the suppression of winter convection (Fig. S1d) associated with427

anomalous surface buoyancy gain. On decadal time scales, sea surface buoyancy (freshwater428

and heat) forcing in the North Atlantic deep convection sites is one of the key factors driving

429

AMOC variability in model simulations (Timmermann et al. 1998; Biastoch et al. 2008; Yeager

430

and Danabasoglu 2014; Pillar et al. 2016). High BBA emissions lead to reduced sea ice

431

formation (Figs. 2d and S1f) and thereby a net buoyancy gain due to fresher surface waters to432

the north of 75°N in the Arctic, while at the same time over 65°N-75°N there is a net buoyancy

433

loss with the decrease in surface freshwater due to the reduced supply of ice (Figs. S3a and434

S3b). However, in the regions where deep convection occurs in this model (the Labrador Sea

435

and the Irminger Sea, Fig. S1c), the buoyancy gains due to surface heating overcome the436

freshwater component of the buoyancy flux anomalies (Fig. S3). As a result, surface buoyancy

437

forcing is positive (buoyancy gain), implying enhanced upper ocean stratification and438

suppressed winter convection that can lead to a slowdown of the AMOC.

439
20
File generated with AMS Word template 2.0

4402) PACIFIC

441

The Pacific STC anomalies (i.e., an anomalous CEC structures) correspond well to zonal442

wind anomalies, as documented in previous studies with idealized atmosphere-ocean coupled443

models (Green and Marshall 2017). During the period when CMIP6 and SMBB forcings differ,444

the trade wind anomalies are largest at the latitudes of their climatological maxima

445(approximately 18°N and 18°S), with a weakening in the NH and strengthening in the SH (Fig.

4468a). As a result, meridional Ekman heat transport calculated from the zonal wind stress447

anomalies and ocean temperature (Sato et al. 2002) show southward anomalies in both the NH

448

and SH, except for latitudes close to and on the equator (Fig. 8b). This is consistent with the449

MHT anomalies and the zonal mean circulation anomalies seen in Figs. 7a and 7c, respectively.

450

This seems to provide favorable conditions for the anomalous southward ocean CEHT to451

compensate for the NH decadal warming induced by the BBA emissions difference between

452

CMIP6 and SMBB. Indeed, previous studies (Green and Marshall 2017; Kang et al. 2018) have453

proposed that mechanistic coupling with changes in the trade winds (i.e., Hadley circulation)454

could form an ocean CEC through Ekman dynamics, and create ocean CEHT anomalies in the455same direction as the atmosphere. In contrast, such a relationship is not apparent in the Atlantic

456(Fig. S4).

457

21
File generated with AMS Word template 2.0

458
459
460
461
462
463
464
465

Fig. 8. Pacific 100-member ensemble mean climatology (black lines with upper x-axis) and
sub-ensemble mean difference (CMIP6 minus SMBB, red lines with bottom x-axis) in (a) zonal
mean zonal wind stress (𝜏! , positive eastward), (b) zonally integrated meridional Ekman heattransport (EHT, positive northward), and (c) zonal mean surface wind stress curl ((∇ × 𝝉)" ).
Blue-shaded latitudes represent the subtropics in each hemisphere defined from the zonal meansurface wind stress curl. The white dots on the differences (red lines) indicate statisticallysignificant differences at the 95% confidence level based on the bootstrapping method.

466
467

However, the actual current fields that constitute the Pacific STCs are complicated, and the468

STCs and the associated MHTs are not fully explained by the near-surface Ekman currents

469

alone that are directly driven by local trade winds. The factors that play a dominant role in the470

net northward basin-wide (zonally integrated) Pacific MHT reflects a complex interplay

471

between narrow WBCs and broad near-surface currents away from the western boundary (non-472

WBCs), with the relative contributions being latitude-dependent (Figs. 9a and 9b). Over the22
File generated with AMS Word template 2.0473

latitudes of the North Pacific subtropical gyre (north of 20°N), northward transport by the WBC

474(Kuroshio) contributes dominantly to the net northward MHT. In the lower latitudes of the475

North Pacific comprising the tropical circulation system, including the southward-flowing

476

equatorial WBC (Mindanao Current) but excepting the equatorial waveguide region (3°N–3°S),

477

northward Ekman heat transport in the non-WBC region due to prevailing trade winds forms

478

the net northward MHT. In the SH, the northward MHT by the northward WBCs of the tropical479

circulation system [the New Guinea Coastal Undercurrent (NGCU) and North Queensland

480

Current/Hiri Current (NQC/HC), cf. Cravatte et al. 2011] is dominant equatorward of 18°S,

481

while the non-WBC northward heat transport is dominant in the subtropical gyre (at least,

482~34°S as shown in Fig. 9a), resulting in the net northward basin-integrated MHT over all

483

latitudes for the SH.

484

485
486
487
488
489
490
491
492
493
494

Fig. 9. (a and b) 100-member ensemble mean climatology (2000–2020) of and (c and d)
sub-ensemble mean difference (CMIP6-SMBB) in ocean full-depth meridional heat transport(MHT, positive northward), surface wind stress, and zonally integrated ocean MHT. 100ensemble mean climatology of SSH (cm) is shown as contours in (a). The 50–50 difference ofthe meridional Ekman heat transport (EHT, same as Fig. 8b) are also shown as the red line in(d). In (c), statistically insignificant MHT differences at the 95% confidence level are crosshatched and only wind stress differences for which the zonal components are statisticallysignificant at the 95% confidence level are shown. Statistically significant zonally integrated
MHT differences are shown as thick solid lines in (d).

49523File generated with AMS Word template 2.0

496

As was the case with the climatological currents, the southward MHT anomalies due to the497

BBA-induced anomalous Pacific CEC reveal latitude-dependence (Figs. 9c and 9d). In the

498

North Pacific subtropical gyre (north of 20°N), anomalous ocean meridional heat flux (MHF)499

by the WBCs is the principal contribution to the Pacific basin’s net southward MHT anomalies

500(Figs. 10a). The northward MHT by the North Pacific subtropical WBC (the Kuroshio) is

501

weakened, such that the zonally-integrated MHT anomalies are southward over the respective

502

latitudes. On the other hand, over the equatorward portion of the South Pacific subtropical gyre503

(18°S–30°S), anomalous southward MHF due to a strengthened South Pacific subtropical504

WBC (the Eastern Australian Current) contributes to the basin-wide southward MHT

505

anomalies over the respective latitudes (Figs. 10e). These WBC transport changes are driven

506

by anomalies in the surface wind stress curl and its zonal integral (Fig. 8c). In other words, the507North Pacific subtropical gyre weakens and the Southern Pacific subtropical gyre strengthens.

508

24File generated with AMS Word template 2.0

509
510
511
512
513
514

Fig. 10. The difference in CMIP6 and SMBB ensemble means of meridional ocean heat
flux (MHF, positive northward, shadings) along five latitudes [a: 20°N, b: 10°N, c: Equator(2°S–2°N), d: 10°S, and e: 20°S]. Upper panels on each depth-longitude section show eastward
cumulative sum of depth-integrated MHF 50–50 difference from the western boundary (i.e.,their values at the eastern boundary are identical to Fig. 9d). Climatological mean (10025
File generated with AMS Word template 2.0

515
516
517
518
519

ensemble mean) MHFs (106 W m-2) are superimposed in the lower panels as contours [contour
levels = (±100., ±10., ±1., ±0.5, 0.)]. Statistically insignificant MHF differences at the 95%
confidence level are cross-hatched. The green lines in the (a), (b), (d), and (e) indicate
approximate Ekman layer depth defined as '2𝐴# /𝑓, where 𝐴# and 𝑓 are the eddy viscosity(𝐴# = 0.05 m2 s-1) and the Coriolis parameter, respectively.

520
521

In contrast to the subtropics, over the latitudes associated with the tropical circulation522

system, MHT anomalies promoted through the Ekman heat transport in non-WBC regions play523

an important role in generating the net southward Pacific basin MHT anomalies (Fig. 9d). In

524

the North Pacific component of the tropical Pacific circulation (2°–18°N), reduced northward

525

Ekman heat transport due to weakened trade winds to the east of 130°E contributes to the net

526

southward zonally-integrated MHT anomaly (Fig. 10b). For the South Pacific component of527

the tropical Pacific circulation (2°–16°S), the net southward zonally-integrated MHT

528

anomalies are caused by enhanced southward Ekman heat transport due to trade wind529

strengthening, a characteristic that is widely distributed over the ocean basin at 10°S (Fig. 10d).

530

At the equator, the response is more complicated. MHT anomalies associated with the531

enhanced narrow southward current around 155°E mainly contribute to the net southward

532

CEHT anomaly, with secondary contributions from near-surface southward CEHT anomalies

533

east of 130°W (Fig. 10c). This equatorial behavior, which demonstrates a dominant

534

contribution of the frictional WBC to sustain anomalous CEHT at the equator, is more535

consistent with recent results obtained with a low-resolution idealized atmosphere-ocean536

coupled model by Green et al. (2019) than it is with the framework whereby the CEHT is537

sustained by a shifting of the STCs themselves (Schneider 2017). That is to say, from the538

northern boundary of the anomalous CEC (around 30°N) to the southern boundary (near 30°S),

539

the dominant contributions for generating the net southward MHT anomalies differs between

540

the contributions in WBC regions and non-WBC regions, depending on their latitude.

541

The strong convergence of the southward zonally-integrated MHT anomalies south of 10°S

542

(Fig. 9d) implies a downward shunting of heat into the ocean interior, thereby accounting for543

the accumulation of heat in the low-latitude Pacific (Equatorial and South Pacific) already544

identified in Fig. 6g. For the South Pacific subtropics (from approximately 18°S to 40°S), with545

the exception of the coastal regions, vertical heat flux (VHF) just below the layer where it

546interacts directly with the atmosphere (below about 100 m depth) corresponds well to what is547

expected from the surface wind stress curl ((∇ × 𝝉)" ) field through Ekman dynamics (Figs. 11a
26
File generated with AMS Word template 2.0

548

and 11b). During the period when CMIP6 and SMBB forcing is different, differences between

549

two 50-member sub-ensembles indicates a positive wind stress curl anomaly over the South

550

Pacific subtropical gyre (Fig. 11c). These positive anomalies enhance the Ekman pumping

551(downwelling) and hence induce anomalous VHFs that carry more heat into the subsurface

552

layer (Fig. 11d). The heat transported below 300 m depth accumulates within the South Pacific553

subtropical gyre, and through time is advected along a circulation structure extending to554

approximately 10°S (Figs. 11e–11g). Part of the heat is then transported to the subsurface zone555

of the equatorial band via the northward NQC/HC and NGCU (Figs. 11h–11j). These heat556

anomalies accumulated within the South Pacific Subtropical gyre and those reached below the557

equatorial thermocline both contribute to the low-latitude Pacific OHC anomalies that persist

558

until the end of the simulation (Figs. 6g and S5). Over a longer-timescale than is the main focus559

of this study, baroclinic wave adjustment associated with an abrupt AMOC change

560(Timmermann et al. 2005) can also be expected to induce subsurface temperature anomalies in561

the tropical and South Pacific (Wang et al. 2018). Here, however, we emphasize the importance562

of the local wind forcing (positive wind stress curl anomalies) in this shorter time-scale ocean563

response.

564

565
566
567
568
569

Fig. 11. 100-member ensemble mean climatology (2000–2020) of (a) surface wind stresscurl ((∇ × 𝝉)" ) and SSH [contours, (cm)] and (b) 100-300 m mean vertical heat flux (VHF,
positive downward) and streamlines at 300 m depth. Differences in sub-ensemble means(CMIP6 minus SMBB) in (c) (∇ × 𝝉)" , (d) VHF, and (e–j) subsurface ocean heat content (OHC,
27
File generated with AMS Word template 2.0

570
571

300–1000 m). Statistically insignificant differences at the 95% confidence level are crosshatched in (c–j).

572
573

5. Summary and Discussion574

We have investigated the ocean response to anomalous NH warming induced by a rectified

575

response to interannual modulations of BBA emissions and corresponding shortwave radiation.

576

We compared two 50-member sub-ensembles of the 100-member CESM2-LE, with their BBA

577

emission forcing differing over 1990–2020. By virtue of having two large sets of simulation

578

with 50 ensemble members in each, which were conducted by a fully coupled Earth system579

model, our experimental design facilitates a clear separation of the forced response to changes

580

in the interannual variability of BBA emissions from background natural variability. The main581

takeaways of this study are (1) the detailed description of the complexity of the Pacific CEC

582

response and (2) the explanation of the mechanism for heat accumulation at the Pacific583

subsurface, resulting in the long-lasting oceanic memory of the interannual BBA forcing. Our

584

emphasis on the ocean response to hemispherically asymmetric thermal forcing, especially in585

the anomalous Pacific CEC, is complementary to previous more atmospherically-focused586

studies that have applied atmosphere-slab ocean coupled models (e.g., Manabe and Broccoli

587

1985; Chiang et al. 2003; Kang et al. 2008), and studies with idealized atmosphere-ocean588

coupled models (e.g., Green and Marshall 2017) , as well as studies applying state-of-the-art

589

climate models with limited numbers of ensemble members (e.g., White et al. 2018; Yu and590

Pritchard 2019), along with studies that have considered longer-timescale responses (from591interdecadal to equilibrium) responses (e.g., Tomas et al. 2016; Wang et al. 2018, 2022).

592

Detectable differences in the responses of the Pacific and Atlantic Oceans are summarized593

schematically in Fig. 12. The initial atmospheric and sea-ice responses to anomalous heat input

594

through radiation across the top of the atmosphere in the NH high latitudes are common to both595

basins (e.g., atmospheric warming, increased in atmospheric southward CEHT, and sea ice

596

melting), but the subsequent ocean responses are notably different between the two basins. In

597

the Atlantic, heat absorbed by the ocean warms the upper part of the North Atlantic and598

enhances upper ocean stratification (Fig. S6a). Some of the heat is subducted into the599

subsurface within the subtropical and subpolar gyres (Fig. S6b) and the enhanced stratification

600

slows the AMOC down. The warm anomalies in the upper North Atlantic Ocean are

601

compensated by the reduction in northward ocean heat transport due to the weakened AMOC
28
File generated with AMS Word template 2.0602

which occurs broadly across the basin. Subsequent to the time interval when forcing is603

divergent, part of the forced warm anomaly remains in the subsurface for multiple decades

604(Figs. S6c–S6f), while the ocean continually releases anomalous heat into the atmosphere.

605

606
607
608
609
610
611
612
613

Fig. 12. Schematic figure for illustrating (a) Pacific and (b) Atlantic Ocean responses to the
NH decadal warming induced by the high-amplitude BBA emissions variability. The top panels
show differences in sub-ensemble mean (CMIP6 minus SMBB for 2000–2020) net surface heat
flux (positive downward) and surface wind as shadings and vectors, respectively. The contours
in the top panel indicate 100-member mean SSH. In the side panels, differences in subensemble mean temperature (shadings) and MOC (contours) are shown. Note that only
negative MOC anomalies from Fig.7 are shown as contours.

614
615

Like the Atlantic, the upper Pacific Ocean is warmed directly during the period when BBA

616

emissions forcing is different (Fig. S5a). In contrast to the Atlantic, however, the anomalous

617

CEC in the Pacific that is induced is associated with hemispherically asymmetric trade wind618

anomalies, characterized by weakened NH trade winds and strengthened SH trade winds. We

619

found that the anomalous CEC and the resultant net southward ocean MHT anomalies covering

620

the tropical Pacific (30°S–30°N) reflect a complex latitude-dependent interplay between the621

WBCs and the near-surface Ekman transport away from the western boundary regions. The

622

deceleration of the North Pacific subtropical WBC (Kuroshio), a component of the anomalous

623

CEC, reduces northward heat transport into the North Pacific and thereby to an extent partly

624

compensates for the anomalous heat flux into the North Pacific. Other components of the heat

625

flux perturbation applied across the surface of the North Pacific are subsequently released to626

the atmosphere and/or transferred into the subsurface (Figs. S5a–S5c). In the South Pacific, the29File generated with AMS Word template 2.0

627

southward heat transport anomalies converge and by this means heat is shunted into the628

subsurface. Heat accumulated in the subsurface of the South Pacific subtropical circulation is629

then partly transported to the layer below the equatorial thermocline via the NGCU, which

630

contributes to the persistent ocean memory (Figs. S5d–S5f).

631

The primary factors distinguishing the ocean responses between the Pacific and Atlantic

632

are the presence of a deep MOC structure for the Atlantic relative to the Pacific, and the633

enhanced efficiency for the Pacific of zonally-integrated wind stress curl-driven transport

634

anomalies, as a consequence of the significantly broader zonal extent. For the present case, the635

anomalous increase in OHC is smaller in the Atlantic than in the Pacific, even though the636

Atlantic gains 30% more heat through local heat flux anomalies during the period of the high-637

amplitude BBA emissions variability (Fig. S7). For the Atlantic, the AMOC is able to638

compensate NH warming anomalies by virtue of not only its access to a deep heat reservoir but

639

also the broad meridional reach of its coherent circulation pattern, meaning that there is more

640

direct leverage for an anomalous OHTC reduction through its southern end (Fig. S7). As the641

Pacific is a much wider basin than the Atlantic, net CEC anomalies can be more effectively

642

controlled by anomalies in trade wind strength associated with hemispherically asymmetric643

thermal forcing perturbations.

644

This study has broader implications for understanding the response of the climate system645

to BBA emissions. As has been described above, previous studies (Fasullo et al. 2022;

646

DeRepentigny et al. 2022; Kim et al. under review) have shown that the BBA emissions forcing

647

protocols for CMIP6 (BB4CMIP6), with spurious discontinuities in the variance of BBA

648

emissions, are somewhat of a mixed blessing. Despite the increased realism in this forcing

649

during data-rich periods, spurious changes in the amplitude of BBA emissions can lead to650

spurious mean state changes in a number of important climate variables. Viewing this in a651

positive light, this sensitivity to large changes in the variance of BBA emissions has also drawn

652

our attention to a previously under-appreciated role for BBA emissions in the climate system.

653

Moving forward, there are two points that we wish to emphasize that are informed by what654

we have learned thus far. First, the rectified heat anomalies in the atmosphere associated with655

changes in the interannual variability of BBA emissions are only an order of 0.4 W m-2, and656

yet this is sufficient to trigger large-scale adjustment processes reaching into the SH. It is657

valuable to consider this within the context of volcanic eruptions (Church et al. 2005;

658

Stenchikov et al. 2009) where the short-term intensity of forcing from an individual volcano
30
File generated with AMS Word template 2.0

659

can be ten times larger in magnitude (~ -4 W m-2 vs. ~0.4 W m-2 of TOA radiation anomalies),

660

but where the net the Earth system perturbation can be of similar amplitude [maximum global

661

OHC anomaly of ~50 ZJ for a Pinatubo-level volcanic eruption (Stenchikov et al. 2009)

662

compared to ~20 ZJ for increased interannual BBA variability in this study]. Although663

individual TOA anomalies associated with high-amplitude modulations of BBA emissions

664

variability are smaller, the cumulative impact of intermittent anomalies can have impacts on665

the system comparable to those of a large volcanic eruption event.

666

Second, the fact that the forcing perturbation of CMIP6-SMBB over 1990–2020 spans only

667

12% of the full 251 years of the CESM2-LE (1850-2100) implies that if realistic BBA

668

emissions were used over the full 1850-2100 period, one might expect that the cumulative669

ocean heat response to the forcing could be larger than it is under SMBB forcing perturbation670

considered here. This could be particularly pertinent if interannual variability in BBA

671

emissions were to increase in time with climate change (e.g., van der Werf et al. 2017). With

672

this in mind, it is recommended that the work presented here be extended by considering the673

climate feedbacks associated with realistic modulations of the frequency and/or intensity of674

fire occurrence over scenario runs spanning 1850-2100. Such an effort would benefit not only

675

our understanding of the processes that contribute to climate change, but also the development676

of new protocols that can be applied to the CMIP7 modeling project.

677
678

Acknowledgments.

679

This study was supported by the Institute for Basic Science (IBS), Republic of Korea, under

680

IBS-R028-D1. This material is based upon work supported by the National Center for

681

Atmospheric Research, which is a major facility sponsored by the NSF under Cooperative

682

Agreement 1852977. MFS was supported by NSF grant AGS-2141728 and NOAA’s Climate

683

Program Office’s Modeling, Analysis, Predictions, and Projections (MAPP) program grant684

NA20OAR4310445. This is IPRC publication X and SOEST contribution Y. The CESM2

685

simulations were conducted on the IBS/ICCP supercomputer ‘Aleph’, 1.43 peta flops high-686

performance Cray XC50-LC Skylake computing system with 18,720 processor cores, 9.59 PB

687

storage and 43 PB tape archive space. We also acknowledge the support of KREONET.

688
689

Data Availability Statement.

31
File generated with AMS Word template 2.0

690
691

The CESM2-LE data is available through NCAR’s Climate Data Gateway athttps://www.earthsystemgrid.org/dataset/ucar.cgd.cesm2le.output.html.

692
693
694

REFERENCES
Archibald, S., C. E. R. Lehmann, J. L. Gómez-Dans, and R. A. Bradstock, 2013: Defining695

pyromes and global syndromes of fire regimes. Proceedings of the National Academy of696

Sciences, 110, 6442–6447, https://doi.org/10.1073/pnas.1211466110.

697

Biastoch, A., C. W. Böning, J. Getzlaff, J.-M. Molines, and G. Madec, 2008: Causes of698

Interannual–Decadal Variability in the Meridional Overturning Circulation of the699

Midlatitude North Atlantic Ocean. Journal of Climate, 21, 6599–6615,

700

https://doi.org/10.1175/2008JCLI2404.1.

701

Bogenschutz, P. A., A. Gettelman, C. Hannay, V. E. Larson, R. B. Neale, C. Craig, and C.-C.

702

Chen, 2018: The path to CAM6: coupled simulations with CAM5.4 and CAM5.5.

703

Geoscientific Model Development, 11, 235–255, https://doi.org/10.5194/gmd-11-235-

704

2018.

705

Broccoli, A. J., K. A. Dahl, and R. J. Stouffer, 2006: Response of the ITCZ to Northern

706

Hemisphere cooling. Geophysical Research Letters, 33,707

https://doi.org/10.1029/2005GL024546.

708

Chiang, J. C. H., and A. R. Friedman, 2012: Extratropical Cooling, Interhemispheric Thermal

709

Gradients, and Tropical Climate Change. Annu. Rev. Earth Planet. Sci., 40, 383–412,

710

https://doi.org/10.1146/annurev-earth-042711-105545.

711

Chiang, J. C. H., M. Biasutti, and D. S. Battisti, 2003: Sensitivity of the Atlantic Intertropical

712

Convergence Zone to Last Glacial Maximum boundary conditions. Paleoceanography,

713

18, https://doi.org/10.1029/2003PA000916.

714

Church, J. A., N. J. White, and J. M. Arblaster, 2005: Significant decadal-scale impact of715

volcanic eruptions on sea level and ocean heat content. Nature, 438, 74–77,716

https://doi.org/10.1038/nature04237.

717

Cravatte, S., A. Ganachaud, Q.-P. Duong, W. S. Kessler, G. Eldin, and P. Dutrieux, 2011:

718

Observed circulation in the Solomon Sea from SADCP data. Progress in Oceanography,

719

88, 116–130, https://doi.org/10.1016/j.pocean.2010.12.015.
32
File generated with AMS Word template 2.0720Danabasoglu, G., S. C. Bates, B. P. Briegleb, S. R. Jayne, M. Jochum, W. G. Large, S.

721

Peacock, and S. G. Yeager, 2012: The CCSM4 Ocean Component. J. Climate, 25, 1361–

722

1389, https://doi.org/10.1175/JCLI-D-11-00091.1.

723

Danabasoglu, G., and Coauthors, 2020: The Community Earth System Model Version 2

724(CESM2). Journal of Advances in Modeling Earth Systems, 12, e2019MS001916,725

https://doi.org/10.1029/2019MS001916.

726

DeRepentigny, P., and Coauthors, 2022: Enhanced simulated early 21st century Arctic sea ice

727

loss due to CMIP6 biomass burning emissions. Science Advances, 8, eabo2405,

728

https://doi.org/10.1126/sciadv.abo2405.

729

Deser, C., R. A. Tomas, and L. Sun, 2015: The Role of Ocean–Atmosphere Coupling in the730

Zonal-Mean Atmospheric Response to Arctic Sea Ice Loss. Journal of Climate, 28, 2168–

731

2186, https://doi.org/10.1175/JCLI-D-14-00325.1.

732

Eyring, V., S. Bony, G. A. Meehl, C. A. Senior, B. Stevens, R. J. Stouffer, and K. E. Taylor,733

2016: Overview of the Coupled Model Intercomparison Project Phase 6 (CMIP6)734

experimental design and organization. Geosci. Model Dev., 9, 1937–1958,735

https://doi.org/10.5194/gmd-9-1937-2016.

736

Fasullo, J. T., J.-F. Lamarque, C. Hannay, N. Rosenbloom, S. Tilmes, P. DeRepentigny, A.

737

Jahn, and C. Deser, 2022: Spurious Late Historical-Era Warming in CESM2 Driven by

738

Prescribed Biomass Burning Emissions. Geophysical Research Letters, 49,739

e2021GL097420, https://doi.org/10.1029/2021GL097420.

740

Green, B., and J. Marshall, 2017: Coupling of Trade Winds with Ocean Circulation Damps

741

ITCZ Shifts. Journal of Climate, 30, 4395–4411, https://doi.org/10.1175/JCLI-D-16-

742

0818.1.

743
744
745

——, ——, and J.-M. Campin, 2019: The ‘sticky’ ITCZ: ocean-moderated ITCZ shifts. Clim
Dyn, 53, 1–19, https://doi.org/10.1007/s00382-019-04623-5.
Hawcroft, M., J. M. Haywood, M. Collins, and A. Jones, 2018: The contrasting climate746

response to tropical and extratropical energy perturbations. Clim Dyn, 51, 3231–3249,

747

https://doi.org/10.1007/s00382-018-4076-8.

748

Heyblom, K. B., H. A. Singh, P. J. Rasch, and P. DeRepentigny, 2022: Increased Variability

749

of Biomass Burning Emissions in CMIP6 Amplifies Hydrologic Cycle in the CESM2
33
File generated with AMS Word template 2.0

750Large Ensemble. Geophysical Research Letters, 49, e2021GL096868,

751

https://doi.org/10.1029/2021GL096868.

752

Hunke, E. C., W. H. Lipscomb, A. K. Turner, N. Jeffery, and S. Elliott, 2015: CICE: the Los

753

Alamos Sea Ice Model Documentation and Software User’s Manual Version 5. 116,754

http://www.ccpo.odu.edu/~klinck/Reprints/PDF/cicedoc2015.pdf.

755
756
757

Kang, S. M., 2020: Extratropical Influence on the Tropical Rainfall Distribution. Curr Clim
Change Rep, 6, 24–36, https://doi.org/10.1007/s40641-020-00154-y.
——, I. M. Held, D. M. W. Frierson, and M. Zhao, 2008: The Response of the ITCZ to758

Extratropical Thermal Forcing: Idealized Slab-Ocean Experiments with a GCM. Journal759

of Climate, 21, 3521–3532, https://doi.org/10.1175/2007JCLI2146.1.

760——, Y. Shin, and S.-P. Xie, 2018: Extratropical forcing and tropical rainfall distribution:

761

energetics framework and ocean Ekman advection. npj Clim Atmos Sci, 1, 1–10,

762

https://doi.org/10.1038/s41612-017-0004-6.

763——, S.-P. Xie, C. Deser, and B. Xiang, 2021: Zonal mean and shift modes of historical764

climate response to evolving aerosol distribution. Science Bulletin, 66, 2405–2411,

765

https://doi.org/10.1016/j.scib.2021.07.013.

766

Lamarque, J.-F., and Coauthors, 2010: Historical (1850–2000) gridded anthropogenic and767

biomass burning emissions of reactive gases and aerosols: methodology and application.

768

Atmospheric Chemistry and Physics, 10, 7017–7039, https://doi.org/10.5194/acp-10-769

7017-2010.

770

Lawrence, D. M., and Coauthors, 2019: The Community Land Model Version 5: Description

771

of New Features, Benchmarking, and Impact of Forcing Uncertainty. Journal of Advances

772

in Modeling Earth Systems, 11, 4245–4287, https://doi.org/10.1029/2018MS001583.

773

Manabe, S., and A. J. Broccoli, 1985: The influence of continental ice sheets on the climate774

of an ice age. Journal of Geophysical Research: Atmospheres, 90, 2167–2190,

775

https://doi.org/10.1029/JD090iD01p02167.

776

van Marle, M. J. E., and Coauthors, 2017: Historic global biomass burning emissions for

777

CMIP6 (BB4CMIP) based on merging satellite observations with proxies and fire models

778(1750–2015). Geosci. Model Dev., 10, 3329–3357, https://doi.org/10.5194/gmd-10-3329-

779

2017.
34
File generated with AMS Word template 2.0

780

McCreary, J. P., and P. Lu, 1994: Interaction between the Subtropical and Equatorial Ocean

781

Circulations: The Subtropical Cell. Journal of Physical Oceanography, 24, 466–497,

782

https://doi.org/10.1175/1520-0485(1994)024<0466:IBTSAE>2.0.CO;2.

783

Pillar, H. R., P. Heimbach, H. L. Johnson, and D. P. Marshall, 2016: Dynamical Attribution

784

of Recent Variability in Atlantic Overturning. Journal of Climate, 29, 3339–3352,

785

https://doi.org/10.1175/JCLI-D-15-0727.1.

786

Rabin, S. S., and Coauthors, 2017: The Fire Modeling Intercomparison Project (FireMIP),

787

phase 1: experimental and analytical protocols with detailed model descriptions. Geosci.

788

Model Dev., 10, 1175–1197, https://doi.org/10.5194/gmd-10-1175-2017.

789

Rodgers, K. B., and Coauthors, 2021: Ubiquity of human-induced changes in climate

790

variability. Earth System Dynamics, 12, 1393–1411, https://doi.org/10.5194/esd-12-1393-

791

2021.

792

Sato, O. T., P. S. Polito, and W. Timothy Liu, 2002: Intra-decadal variability in the Ekman

793

heat flux from scatterometer winds: INTRA-DECADAL VARIABILITY IN THE

794

EKMAN HEAT FLUX. Geophys. Res. Lett., 29, 12-1-12–14,795

https://doi.org/10.1029/2002GL014775.

796

Schneider, T., 2017: Feedback of Atmosphere-Ocean Coupling on Shifts of the Intertropical

797

Convergence Zone. Geophysical Research Letters, 44, 11,644-11,653,

798

https://doi.org/10.1002/2017GL075817.

799——, T. Bischoff, and G. H. Haug, 2014: Migrations and dynamics of the intertropical

800

convergence zone. Nature, 513, 45–53, https://doi.org/10.1038/nature13636.

801

Smith, R., and Coauthors, 2010: The Parallel Ocean Program (POP) reference manual: Ocean

802

component of the Community Climate System Model (CCSM).

803

http://n2t.net/ark:/85065/d70g3j4h.

804

Stenchikov, G., T. L. Delworth, V. Ramaswamy, R. J. Stouffer, A. Wittenberg, and F. Zeng,805

2009: Volcanic signals in oceans. Journal of Geophysical Research: Atmospheres, 114,806

https://doi.org/10.1029/2008JD011673.

807

Timmermann, A., M. Latif, R. Voss, and A. Grötzner, 1998: Northern Hemispheric

808

Interdecadal Variability: A Coupled Air–Sea Mode. Journal of Climate, 11, 1906–1931,809

https://doi.org/10.1175/1520-0442(1998)011<1906:NHIVAC>2.0.CO;2.
35
File generated with AMS Word template 2.0

810——, S.-I. An, U. Krebs, and H. Goosse, 2005: ENSO Suppression due to Weakening of the811

North Atlantic Thermohaline Circulation. Journal of Climate, 18, 3122–3139,

812

https://doi.org/10.1175/JCLI3495.1.

813

Tomas, R. A., C. Deser, and L. Sun, 2016: The Role of Ocean Heat Transport in the Global

814

Climate Response to Projected Arctic Sea Ice Loss. Journal of Climate, 29, 6841–6859,

815

https://doi.org/10.1175/JCLI-D-15-0651.1.

816

Wang, K., C. Deser, L. Sun, and R. A. Tomas, 2018: Fast Response of the Tropics to an

817

Abrupt Loss of Arctic Sea Ice via Ocean Dynamics. Geophysical Research Letters, 45,818

4264–4272, https://doi.org/10.1029/2018GL077325.

819——, L. Wu, H. Liu, B. Dan, H. Dai, and C. Deser, 2022: Contrary Responses of the Gulf

820

Stream and the Kuroshio to Arctic Sea Ice Loss. Atmosphere, 13, 514,821

https://doi.org/10.3390/atmos13040514.

822

van der Werf, G. R., and Coauthors, 2017: Global fire emissions estimates during 1997–

823

2016. Earth System Science Data, 9, 697–720, https://doi.org/10.5194/essd-9-697-2017.

824

White, R. H., A. A. McFarlane, D. M. W. Frierson, S. M. Kang, Y. Shin, and M. Friedman,825

2018: Tropical Precipitation and Cross-Equatorial Heat Transport in Response to826

Localized Heating: Basin and Hemisphere Dependence. Geophysical Research Letters,

827

45, 11,949-11,958, https://doi.org/10.1029/2018GL078781.

828

Yeager, S., and G. Danabasoglu, 2014: The Origins of Late-Twentieth-Century Variations in829

the Large-Scale North Atlantic Circulation. Journal of Climate, 27, 3222–3247,

830

https://doi.org/10.1175/JCLI-D-13-00125.1.

831

Yoshimori, M., A. Abe-Ouchi, H. Tatebe, T. Nozawa, and A. Oka, 2018: The Importance of

832

Ocean Dynamical Feedback for Understanding the Impact of Mid–High-Latitude

833

Warming on Tropical Precipitation Change. Journal of Climate, 31, 2417–2434,

834

https://doi.org/10.1175/JCLI-D-17-0402.1.

835

Yu, S., and M. S. Pritchard, 2019: A Strong Role for the AMOC in Partitioning Global

836

Energy Transport and Shifting ITCZ Position in Response to Latitudinally Discrete Solar

837

Forcing in CESM1.2. Journal of Climate, 32, 2207–2226, https://doi.org/10.1175/JCLI-

838

D-18-0360.1.

839
36File generated with AMS Word template 2.0