Generated using the official AMS LATEX template v6.1

1

The CESM2 Single Forcing Large Ensemble and Comparison to CESM1:

2

Implications for Experimental Design3Isla R. Simpsona , Nan Rosenblooma , Gokhan Danabasoglua , Clara Desera , Stephen G. Yeagera ,4

Christina S. McCluskeya , Ryohei Yamaguchib , Jean-Francois Lamarquea , Simone Tilmesc ,5Michael J. Millsc , Keith B. Rodgersd,e6

a Climate and Global Dynamics Laboratory, National Center for Atmospheric Research, Boulder,7

Colorado

8b Japan Agency for Marine-Earth Science and Technology, Yokosuka, Japan9c Atmospheric Chemistry Observations and Modelling Laboratory, National Center for

10

Atmospheric Research, Boulder, Colorado

11d Center for Climate Physics, Insitute for Basic Science, Busan, South Korea

12

e Pusan National University, Busan, South Korea

13

Corresponding author: Isla Simpson, islas@ucar.edu
1

14

ABSTRACT: Single Forcing Large Ensembles are a relatively new tool for quantifying the con-15

tributions of different anthropogenic and natural forcings to the historical and future projected16

evolution of the climate system. This study introduces a new single forcing large ensemble with17

the Community Earth system Model version 2 which can be used to separate the influences of18

greenhouse gases, anthropogenic aerosols, biomass burning aerosols, and all remaining forcings,

19

on the evolution of the Earth System from 1850 to 2050. Here, the forced responses of global near

20

surface temperature and associated drivers are examined in CESM2 and compared with those in

21

a single forcing large ensemble with CESM2’s predecessor, CESM1. The experimental design,22

the imposed forcing and the model physics all differ between the CESM1 and CESM2 ensembles.

23

In CESM1 an “all-but-one” approach was used where everything except the forcing of interest is

24

time evolving, while in CESM2 an “only” approach is used, where only the forcing of interest is

25

time evolving. This experimental design choice is shown to matter considerably for anthropogenic

26

aerosol-forced change in CESM2, due to state dependence of cryospheric albedo feedbacks and27

non-linearity in the Atlantic Meridional Overturning Circulation (AMOC) response to forcing.

28

This impact of experimental design is, however, strongly dependent on the model physics and/or

29

the imposed forcing as the same sensitivity to experimental design is not found in CESM1, which

30

appears to be an inherently less non-linear model in both its AMOC behavior and cryospheric31

feedbacks.

232

1. Introduction

33

Historically, the Earth’s climate system has evolved under a mixture of natural and anthropogenic

34

forcings and it will continue to do so moving forward. A common approach that is used to35

disentangle and understand the relative contributions of such forcings to the evolution of the

36

climate system is to perform Earth System Model (ESM) experiments in which only some forcings

37

are evolving in time while others are held fixed. These experiments, which we will refer to as single

38

forcing experiments even though they may be used to isolate the influence of multiple forcings at

39

once, are most informative when a relatively large number of ensemble members are available, such

40

that the forced signal can be isolated from the internal variability (e.g., Deser et al. 2020a). Many

41

modelling centers have performed single forcing experiments under the coordinated framework of

42

the “Detection and Attribution Model Intercomparison Project” (DAMIP, Gillett et al. 2016) as

43part of the Coupled Model Intercomparison Project (CMIP) and their utility has been recognized

44

by the World Climate Research Program through emphasis on single forcing large ensembles as

45part of the lighthouse activity on “Explaining and Predicting Earth System Change” (Smith et al.

462022).

47

Single forcing experiments have been a core component of Intergovernmental Panel on Climate

48

Change (IPCC) reports and have been used to conclude that human influence has unequivocally49

warmed the climate (IPCC 2021; Gillett et al. 2021). Beyond this, they have been used to50

investigate the wide ranging impacts of individual forcings on various aspects of the climate system.

51

These include: the global patterns of surface temperature and precipitation anomalies induced by52

greenhouse gas versus aerosol forcing (Deser et al. 2020b; Shi et al. 2022); the influence of aerosol

53

forcing on precipitation in the Sahel region (Dong et al. 2014; Giannini and Kaplan 2019; Hirasawa

54

et al. 2020; Zhang et al. 2021) and other monsoon regions of the world (Li et al. 2018; Undorf

55

et al. 2018; Monerie et al. 2022); the influence of ozone depleting substances and greenhouse56

gases on precipitation over Australia (Delworth and Zeng 2014); assessment of the aerosol-forced

57

contribution to trends in the Pacific Ocean (Allen et al. 2014; Dittus et al. 2021) ; assessment of the58

counteracting influence of greenhouse gases and aerosols on Arctic sea ice (Mueller et al. 2018)59

and Arctic temperatures (England et al. 2021); the impacts of individual forcings on the North

60

Atlantic Ocean circulation (Watanabe and Tatebe 2019; Dagan et al. 2020; Baek et al. 2022) and61

global sea level rise (Fasullo et al. 2020); the impacts of land cover change and irrigation on surface3

62

temperature and precipitation (Singh et al. 2018); and, the impact of individual forcings on various

63

other aspects of the hydrological cycle and extreme weather (Chiang et al. 2021; Pendergrass et al.

64

2019; Bonfils et al. 2020; Seong et al. 2021; Touma et al. 2021). Single forcing experiments can

65

also provide a useful testbed for exploring model sensitivity to differences in imposed forcings

66(e.g., Fyfe et al. 2021) or for inter-comparing the response to forcings among models (e.g., Menary

67

et al. 2020; Dittus et al. 2021).

68

In the design of single forcing experiments choices must be made. Under the DAMIP protocol,

69

the forcing of interest is evolving in time while all others are held fixed at pre-industrial values,

70

referred to as the “only” method, hereafter. Another option is to evolve all forcings in time except

71

the one of interest and determine that forcings influence by differencing this experiment from an

72

all forcing simulation, referred to as the “all-but-one” method, herafter; this was the choice made73

for the single forcing large ensemble with the CESM1 model. There is also a choice as to what74

year forcings are held fixed at; DAMIP fixes them at 1850, while the CESM1 single forcing large

75

ensemble fixed them at 1920. Whether these various design choices will produce the same answer

76

as to a forcing’s influence will depend on whether there are substantial non-linearities or state77

dependencies in the system, and prior results have been mixed as to whether this is the case. Meehl

78

et al. (2004) found the global mean temperature response to forcings was approximately linearly

79

additive while Feichter et al. (2004) and Ming and Ramaswamy (2009) found that it was not. A

80

more recent study by Deng et al. (2020) assessed additivity of the response of greenhouse gases

81

and aerosols in time slice experiments with CESM1. They found that for global mean temperature,82

the influence of these forcings were approximately linearly additive but for other features, such as

83

autumn Arctic sea ice cover and East Asian precipitation, non-linearities did exist.

84

Here, we present a new single forcing large ensemble with the Community Earth System Model85

version 2 (CESM2). The aims of this study are two-fold: (1) to introduce this new dataset

86

that researchers can use to further probe the impacts of individual forcings on the evolution of87

the Earth system according to this model and (2) to understand differences in the global mean

88

temperature and radiative responses between this single forcing large ensemble and its predecessor

89(the CESM1 single forcing large ensemble (Deser et al. 2020b)). With regards to the second

90

goal, we find substantial differences in the anthropogenic aerosol-forced global mean near surface91

air temperature evolution between the CESM1 and CESM2 ensembles. Three factors have the4

92

potential to contribute to this: differences in imposed aerosol emissions; differences in model

93

physics; and, differences in the experimental design. Here, we make use of additional targeted94

experiments to attempt to isolate the relative role of the experimental design and, while this will95

be shown to have an influence, its impact is found to be sensitive to model physics and/or imposed

96

forcings. The models, experimental design and methods are described in section 2. In section 3,

97

we compare the global mean surface air temperature (GMST) response between the CESM2 and98

CESM1 single forcing large ensembles and reveal a substantial difference in the aerosol-forced

99

response. In section 4, we then explore the influence of experimental design on the aerosol-forced

100

GMST change in the CESM2 single forcing large ensemble and follow this with a comparison to101

CESM1 in section 5. Discussion and conclusions are provided in section 6.

102

2. Models, Experiments and Methods

103

a. CESM2 and its single forcing experiments

1041) The Model105

CESM2 is the latest generation Earth System Model developed by the U.S. National Center106

for Atmospheric Research in collaboration with others (Danabasoglu et al. 2020). The default

107

configuration of CESM2, which was used to contribute experiments to CMIP6 (Eyring et al.

1082016), simulates the global coupled Earth system at approximately 1◦ horizontal resolution. The

109

atmospheric component is the Community Atmosphere Model version 6 (CAM6, Bogenschutz

110

et al. 2018) with a model top at ∼40 km and 32 layers in the vertical. It is coupled to the Parallel

111

Ocean Program version 2 (POP2) ocean model (Smith et al. 2010; Danabasoglu et al. 2012),

112

the Community Land Model version 5 (CLM5, Lawrence et al. 2019) and the Community Ice

113

Code version 5 (CICE5, Hunke et al. 2015) and all the simulations in this study have fixed ice

114

sheets. We refer readers to Danabasoglu et al. (2020) for more details and to the following studies

115

for evaluation of various aspects of CESM2: Lawrence et al. (2019) for the representation of116

land surface processes; Simpson et al. (2020) for the large scale atmospheric circulation and its

117

variability; Meehl et al. (2020) for the representation of monsoons; Capotondi et al. (2020) for118

the representation of Pacific sea surface temperature variability; and, DuVivier et al. (2020) for119

the representation of sea ice. In this description, we focus on the aspects of CESM2 that are of120

particular relevance to the single forcing large ensemble.
5

121

Within CESM2 with CAM6, atmospheric greenhouse gas concentrations are prescribed as

122

monthly time evolving global concentrations. Aerosol forcing is introduced into the model via

123

emissions of black carbon (BC), particulate organic matter (POM), sulfur dioxide (SO2 ), sulfate

124(SO4 ) and secondary organic aerosol precursor gas (SOAG). In CAM6, the aerosol scheme is the125four-mode Modal Aerosol Module (MAM4, Liu et al. 2016). This consists of a very simple sec-

126

ondary organic aerosol scheme that does not include the oxidation of Volatile Organic Compounds

127

and it is not interactively coupled to biogenic emissions (Tilmes et al. 2019). For carbonaceous

128

aerosols, compared to its predecessor (MAM3, Liu et al. 2012), MAM4 contains an additional

129

primary accumulation carbonaceous aerosol mode to allow for an explicit treatment of the micro-

130

physical aging of primary carbonaceous aerosols. Hydrophobic BC and POM are emitted into

131

this fourth primary aerosol mode, where they do not activate cloud condensation nuclei (CCN)132

and cannot be removed by wet deposition. Over time (of the order 2-3 days) they move into the133hydrophilic accumulation mode, where they are available for cloud droplet activation as CCN and

134

also participate in wet deposition. This explicit treatment of the aging of carbonaceous aerosol in

135

MAM4 compared to MAM3, where ageing was instantaneous, has the overall effect of increasing136

the lifetime and subsequent burdens of BC and POM (Liu et al. 2016, and also discussed in the

137

Appendix). The cloud microphysics scheme is version 2 of the Morrison-Gettelman scheme (MG2,138

Gettelman and Morrison 2015) which, unlike its predecesor in CESM1, now includes dependence

139

of mixed-phase immersion freezing ice nucleation on aerosols, i.e., dust aerosol acts as ice nucle-

140

ating particles on which super cooled liquid water or vapor can freeze. CAM6 does not have a141

prognostic representation of stratospheric or tropospheric ozone or volcanic aerosol and, therefore,

142

these forcings are prescribed.

143

In CLM5, each grid cell is composed of multiple land units (vegetated, lake, urban, glacier and144crop) and each land unit has a specified number of columns which are then divided up into multiple

145

patches. These patches contain a plant or crop functional type (PFT or CFT) which is prescribed146

through a land use time series file. Land use and land cover change can, therefore, be introduced

147

through specified evolution of the PFTs and CFTs within each land unit and/or varying the fractional

148

area covered by the land unit components, which allows transitions between natural vegetation,

149

crop and glacier land units (a new feature within CLM5). The land use time series specification also

150

determines the soil texture, wood harvest, industrial Nitrogen fertilizer application amounts and the

6

151

area of the land surface equipped for irrigation. Irrigation is applied dynamically within the model152

to the irrigation-equiped area and is applied to achieve a target soil moisture level (Lombardozzi

153

et al. 2020). When using biogeochemistry mode (as in the experiments here), leaf area index and154

canopy height are prognosed by the model. While CLM5 includes the simulation of fire internally,155

the current default is that this does not produce emissions that are seen by the atmosphere - biomass

156

burning emissions to the atmosphere are prescribed from forcing datasets.

157

2) Forcings and Experiments

158

The CESM simulations used in the following analysis are summarized in Table 1. The baseline159

ensemble for the CESM2 single forcing large ensemble is the second set of 50 members of the160

CESM2 large ensemble (Rodgers et al. 2021). The CESM2 large ensemble, referred to as LENS2

161

hereafter, is a 100-member ensemble of simulations run under CMIP6 historical forcings between

162

1850 and 2014 and forcings of the Shared Socio-economic Pathway 3-7.0 (SSP3-7.0, Meinshausen

163

et al. 2020) thereafter. A mixture of “macro” and “micro” initialization strategies were used

164

to introduce ensemble spread where “macro” refers to initializing each model component from165

different years of the CESM2 pre-industrial control and “micro” refers to introducing ensemble

166

spread through a round-off level perturbation applied to the initial atmospheric potential temperature

167

field. The first and second set of 50 members of LENS2 are run with different biomass burning

168

aerosol emissions over the period 1990 to 2020 (see (Rodgers et al. 2021). The first 50 members

169

of LENS2 use the default CMIP6 biomass burning aerosol dataset, which contains higher levels of

170

interannual variability during 1997-2014 compared to earlier and later periods due to the inclusion

171

of satellite derived emissions. To avoid these artificial discontinuities in biomass burning aerosol

172

variability, which has been shown to cause a rectified climate response in some regions (Fasullo

173

et al. 2021; DeRepentigny et al. 2022), the second set of 50 members of LENS2 were run with a174smoothed (11-year running mean) version of the biomass burning dataset which alters the emissions

175

from 1990 to 2020. All of the CESM2 single forcing large ensemble experiments use the smoothed176

biomass burning emissions dataset, and are compared with the corresponding second 50-member

177

set of LENS2 simulations (see supplemental Fig. 1 for a comparison of the smoothed and default

178

CMIP6 biomass burning emissions).
7

179

180

Table 1. Summary of CESM experiments. The presence of a "1" or "2" in the experiment name indicates thesimulation was performed with CESM1 or CESM2, respectively.

Name

PeriodModel# Members

Description

LENS2

1850-2100

CESM2

50

All forcings evolving. CMIP6 Historical and SSP3-7.0 but with smoothed BMB

AAER2

1850-2050

CESM2

15

Only anthropogenic aerosols evolving. Other forcings fixed at 1850

GHG2

1850-2050

CESM2

15

Only greenhouse gases evolving. Other forcings fixed at 1850

BMB2

1850-2050

CESM2

15

Only biomass burning emisssions evolving. Other forcings fixed at 1850

EE2

1850-2050

CESM2

15

Forcings other than AAER, GHG and BMB evolving. AAER, GHG and BMB fixed at 1850s

LENS1

1920-2100

CESM1

40

All forcings evolving. CMIP5 Historical and RCP8.5 forcings

XAAER1

1920-2080

CESM1

20

All forcings evolving except anthropogenic aerosols which are kept fixed at 1920s levels

XGHG1

1920-2080

CESM1

20

All forcings evolving except greenhouse gases which are kept fixed at 1920s levels

XBMB1

1920-2030

CESM1

15

All forcings evolving except biomass burning which is kept fixed at 1920s levels

XAAER2

1920-2050

CESM2

3

All forcings evolving except anthropogenic aerosols which are kept fixed at 1920s levels

AAER1

1850-2050

CESM1

3

Only anthropogenic aerosols evolving. Other forcings fixed at 1850

181

As summarized in the top portion of Table 1, the CESM2 single forcing large ensemble consists of182

four sub-ensembles of 15 members each that run from 1850 to 2050 following the “only” approach

183

where only the forcing(s) of interest are evolving in time and others are held fixed at 1850s values:184

the greenhouse gas ensemble (GHG2); the anthropogenic aerosol ensemble (AAER2); the biomass

185burning aerosol ensemble (BMB2) and the “Everything Else” ensemble (EE2), where the number

186

2 in each acronym refers to CESM2. The ensemble members differ through a macro initialization,

187

i.e., they are initialized from different years of the pre-industrial control simulation (the same188

years as members 1-10 and 91-95 of LENS2) to minimize any effects of ocean persistence on

189ensemble spread from the beginning of the run. In GHG2, only greenhouse gas concentrations

190

are evolving in time and all other forcings are held fixed at their 1850s values. Fig. 1a shows the191

time-evolution of two of the important greenhouse gases (CO2 and CH4 ) and others can be seen in

192

supplemental Fig. 2. In AAER2, only “anthropogenic aerosol” emissions are evolving in time. By

193“anthropogenic aerosols” here, we refer to industrial, agricultural, domestic and transport related194

emissions and acknowledge that this is not all anthropogenic emissions because it does not include

195

anthropogenic influences on biomass burning (e.g., van Marle et al. 2017). The global emissions of

196

three of the main aerosols or aerosol pre-cursors (SO2 , BC and SO4 ) in the AAER2 ensemble are

197

shown in Fig. 1b and others can be seen in supplemental Fig. 3. In BMB2 only biomass burning

198

emissions are evolving in time, with three of the main emissions sources shown in Fig. 1c and8199

others shown in supplemental Fig. 1. The transitions in variance of biomass burning emissions

200

between the historical portion, the smoothed 1990-2020 period, and then the subsequent SSP3-7.0

201

projection period are quite apparent and this should be improved upon in any emissions datasets

202

that are developed in the future. All other forcings, aside from those that are time evolving in

203

GHG2, AAER2 and BMB2 are time evolving in the EE2 ensemble. Some of the main forcings that204

are evolving in this ensemble are the solar insolation (Fig. 1d black), the stratospheric volcanic

205

aerosol (Fig. 1d red), stratospheric and tropospheric ozone (Fig. 1e for Southern Hemisphere (SH)206

stratospheric ozone) and land use and land cover change (Fig. 1f). The stratospheric volcanic

207

aerosol and tropospheric and stratospheric ozone concentrations are derived from the average of a208

three-member ensemble of simulations with the Whole Atmosphere Community Climate Model

209

version 6 (WACCM6, Gettelman et al. 2019), also run under historical and SSP3-7.0 forcings.

210

This choice was made for all CESM2 experiments to limit forcing differences between CESM2

211

and CESM2-WACCM and, given the similarity between CESM2 and CESM2-WACCM in the212

troposphere and lower stratosphere, the CESM2-WACCM ozone and volcanic aerosol fields will213

be more consistent with the model dynamics and atmospheric structure than the CMIP6 forcing

214

datasets. Overall, each forcing is time evolving in one of these sub-ensembles, allowing the215

additivity of forcing contributions to be tested.

216

The CESM2 pre-industrial control simulation is also used to examine the behavior of the Atlantic

217

Meridional Overturning Circulation (AMOC). This simulation is run for 2000 years under forcings

218

that are representative of 1850s conditions, following the CMIP6 protocol and we make use of219

simulation years 400 to 2000.

227

b. CESM1 and its single forcing experiments

228

1) The Model229

CESM1 is the previous generation of CESM released in 2010 (Hurrell et al. 2013). It has been230

widely used, including through the CESM1 large ensemble (Kay et al. 2015) and the CESM1 single

231

forcing large ensemble (Deser et al. 2020b). Between CESM1 and CESM2, major developments

232

were undertaken in the atmosphere and land components, in particular. The atmospheric component

233

of CESM1 is CAM5 and this contains the older aerosol scheme, MAM3, which, as discussed above,234

does not allow for the explicit treatment of the aging of primary carbonaceous aerosols although
9

400

1.6

350

1.2

300
0.8
1850 1875 1900 1925 1950 1975 2000 2025 2050(c) Biomass burning aerosols (BMB2)
10

100

8

80

6

BC or SO4 (Tg/y)

2.0

450

SO2
BC
SO4

60

4

40

2

20

2.4

SO2
BC
SO4

11

2.2
2.0

10

1.8

9

1.6

0
0
1850 1875 1900 1925 1950 1975 2000 2025 2050

12
SO4 (×10 2 Tg/y)120(b) Anthropogenic aerosols (AAER2)
SO2 or BC (Tg/y)

140

2.4

CH4 (ppmv)

CO2 (ppmv)

500(a) Greenhouse gases (GHG2)CO2CH4SO2 (Tg/y)550

81850 1875 1900 1925 1950 1975 2000 2025 2050

10

340.4

8

340.2
340.0

Solar
Volcanic

64201850 1875 1900 1925 1950 1975 2000 2025 2050(e) 70hPa 60 S-90 S O3 Concentration
2.4
2.3
2.2
2.1
2.0
1.9
1.8
1.7
1.61850 1875 1900 1925 1950 1975 2000 2025 2050

(f) Land cover
Area (million km2)

340.6

Stratospheric AOD (×10 2)O3 (×10 6 mol/mol)

Solar Insolation (Wm 1)

Everything Else (EE2)(d) Solar + Volcanic

TREES
175
150
125
GRASS
100
75
50
SHRUBS
25
CROPS
01850 1875 1900 1925 1950 1975 2000 2025 2050

226

Fig. 1. Annual means of selected forcings and their evolution in the single forcing ensembles. (a) Global CO2 (left axis) andCH4 (right axis) concentrations as they evolve in the GHG2 ensemble. (b) Global emissions of SO2 (left axis) and BC and SO4(right axis) as they evolve in the AAER2 ensemble. (c) Global emissions of SO2 and BC (left axis) and SO4 (right axis) as they
evolve in the BMB2 ensemble. (d) - (f) show various forcings that evolve in the “Everything Else” ensemble: (d) Solar insolation(left axis) and stratospheric Aerosol Optical Depth (AOD) from the WACCM6 simulations that produced the volcanic aerosol
forcing (right axis); (e) 70hPa ozone concentration averaged over 60◦ S to 90◦ S; and (f) surface area covered by various land covertypes (blue hatching shows the irrigated land surface area) and note that bare ground is not shown.

235

the representation of secondary organic aerosols is the same as in MAM4. The cloud microphysics

236

scheme is the first version of the Morrison-Gettelman scheme (MG1, Morrison and Gettelman2372008) which does not relate mixed-phase imersion freezing ice nucleation to aerosols. In CAM6,

238

many other atmospheric parameterizations underwent development compared to CAM5. For

239

example, CAM5 uses an older generation parameterization of shallow convection and boundary

240

layer turbulence compared to CAM6 and it has a simpler representation of orographic drag.

220
221
222
223
224
225

241

The land component of CESM1 is CLM4 (Lawrence et al. 2011). As with the atmosphere, the242

CLM5 model component of CESM2 contains major developments compared to this older model243

version. Major updates were performed on the representation of soil and plant hydrology, snow

244

density, river modeling, carbon and nitrogen cycling and crop modeling. These updates are, in245

general, found to lead to improvements in the representation of many land surface processes in246

the newer generation of the model (Lawrence et al. 2019) with impacts on the representation of247

climate variability (e.g., Simpson et al. 2022). Of relevance for the prescription of time evolving

248

forcings, CLM4 does not allow for the time-evolution of the fractional area of each grid point249

covered by different land cover types, e.g., it does not allow for transitions in the weighting of
10

250

natural vegetation versus crops within a grid cell, just time-evolution of the plant functional types

251

within a given land cover type. It has a much simpler representation of crops and it does not have252

a representation of irrigation.

2532) Forcings and Experiments

254

Forcings within CESM1 are prescribed in a similar manner to CESM2. Greenhouse gases

255

are represented by prescribed global surface concentrations and aerosols are introduced through256

emission sources. Solar variability is introduced through variations in the total solar irradiance,

257

time evolving volcanic aerosols and ozone concentrations are prescribed and land use and land cover258

change is introduced via time evolving plant functional types. CMIP5 era forcings are used in the259

CESM1 experiments. Historical forcings are used prior to 2005 and forcings of the Representative

260

Concentration Pathway 8.5 (RCP8.5) are used thereafter (Meinshausen et al. 2011; Lamarque et al.

2612011). Much like for CESM2, the ozone forcing in the CESM1 experiments is prescribed based262

on WACCM simulations (in this case using WACCM4 and a two member ensemble with a 10-year263

running mean applied to each month of ozone forcing separately). In CESM1 the volcanic aerosol

264

forcing does not come from WACCM, but rather CMIP5 forcing is used.

265

As summarized in the second portion of Table 1, the baseline ensemble for CESM1 is the CESM1

266

large ensemble (Kay et al. 2015), referred to as LENS1 hereafter. This is a 40-member ensemble

267

that runs from 1920 to 2100 in which all forcings are evolving, with each member differing

268

through a round-off level perturbation introduced to the atmospheric potential temperature field269

at initialization (micro initialization). The CESM1 single forcing large ensemble (Deser et al.

270

2020b) consists of 3 sub-ensembles of 15 or 20 members that use the “all-but-one” method where271

all forcings are time evolving except the forcing of interest, which is held fixed at the values for272

1920. Ensemble spread is also introduced through micro initialization and we refer to these sub-273

ensembles as XFORCING where X denotes that FORCING is held fixed. In the XGHG1 ensemble274

all forcings are evolving except greenhouse gases, in the XAAER1 simulation all forcings except

275

anthropogenic aerosols are evolving (again, by anthropogenic aerosols here, we do not include276

anthropogenic influences on biomass burning) and, in the XBMB1 ensemble all forcings except

277

biomass burning are evolving. The number 1 here denotes that these are CESM1 simulations. The

278

time-evolution of the biomass burning emissions, greenhouse gas concentrations and anthropogenic
11

279

aerosol emissions in CESM1 can be compared with those in CESM2 in supplementary Figs. 1 to 3.

280

Note that a fourth ensemble was originally included in which all forcings except land use and land281

cover change were evolving but this dataset has since been retracted due to an error. The XGHG1282

and XAAER1 ensembles extend from 1920 to 2080 and the XBMB1 ensemble extends from 1920283

to 2030. Given that the “all-but-one”approach is used, the influence of a given forcing must be

284

determined by taking the difference between LENS1 and the XFORCING ensemble and because285

not all forcings are represented within an XFORCING ensemble, a complete test of additivity

286

cannot be performed.

287

Years 400 to 2200 of the CESM1 pre-industrial control simulation are also used. This simulation288

was run under forcings that are representative of 1850s conditions, following the CMIP5 protocol.

289

c. Experimental design sensitivity tests

290

As will be shown, substantial differences in the inferred response to aerosol forcing are found291

between CESM1 and CESM2. To test the influence of the method used, i.e., “only” versus

292“all-but-one”, we perform 3-member ensembles of an “only” anthropogenic aerosol experiment293

with CESM1, referred to as AAER1 and an “all-but-one” anthropogenic aerosol experiment with

294

CESM2, referred to as XAAER2, as summarized in the bottom portion of Table 1. In AAER1,

295

CESM1 is run with time evolving anthropogenic aerosol forcing from 1850 to 2050 with all other296

forcings held fixed at 1850s values and members differ through micro initialization. In XAAER2,

297

CESM2 is run from 1920 to 2050 with all forcings evolving except anthropogenic aerosols which298

are held fixed at those of 1920 and the members differ via macro initialization.

299

d. Methods

300

We focus on the period 1920 to 2050, which is common to the majority of CESM simulations

301(Table 1), and, unless otherwise stated, we consider ensemble mean anomalies from the 1920-1940

302

average, which is at the beginning of the CESM1 single forcing large ensemble. For CESM1 the303

influence of a given forcing at a given time period is therefore given by the difference between

304

LENS1−XFORCING at that time period and LENS1−XFORCING averaged over 1920 to 1940.

305

For CESM2, the influence is simply given by the difference between that time period and the306

average of 1920 to 1940 for the single forcing ensemble.
12

307

To quantify uncertainty and statistical significance a bootstrapping approach is used whereby,

308

within each ensemble, members are randomly sampled with replacement and a new ensemble mean309

calculated. This is repeated 1000 times and the uncertainty range on anomalies for that ensemble

310

is given by the 2.5th to 97.5th percentile range of those bootstrapped ensemble means. For the311

3-member sensitivity tests, rather than performing bootstrapping on the 3-member ensembles,312

we either compare them to uncertainty ranges on the larger ensembles that are calculated by313

subsampling 3 members with replacement or we estimate uncertainty and significance levels314

by bootstrapping equivalent sample sizes from the respective pre-industrial control, under the315

assumption that the internal variability of the 1850s climate is representative of that throughout the316

simulation.

317

3. Global mean surface air temperature evolution: CESM2 versus CESM1

318

We begin by comparing the evolution of global mean near surface (2-m) air temperature ([Ts ],319

and we will use [x] throughout to denote the global mean of variable x) between the CESM2 single

320

forcing large ensemble and CESM1 single forcing large ensemble in Fig. 2. Recall that the CESM1

321

single forcing large ensemble does not have all forcing contributions represented, as the equivalent322

of the “Everything Else” ensemble was not performed. We, therefore, estimate the CESM1

323“Everything Else” contribution as the residual LENS1−( (LENS1−XGHG1) + (LENS1−XAAER1)

324+ (LENS1−XBMB1)) for comparison with EE2. This can only be done out to 2030 when the325

XBMB1 simulation ends and assumes linearity, which may not be valid.

334

Figure 2a shows the time-evolution of [Ts ] for LENS2 (black) and the contributions that are

335

inferred to be due to the different forcing components. This can be compared with the equivalent

336

for CESM1 in Fig. 2b. In both CESM1 and CESM2 greenhouse gases (red) act to increase337[Ts ], while anthropogenic aerosols (blue) act to decrease it. The role of biomass burning aerosols

338(brown) in [Ts ] evolution is fairly minimal, while “Everything Else” (green) acts to cool the planet

339

during the major volcanic eruptions of the 20th Century (e.g., El Chichon in the early 1960s and340

Pinatubo in the early 1990s) and to warm the planet relative to 1920-1940 throughout the first half341

of the 21st Century (most apparent in CESM2 where the EE contribution can be examined beyond

3422030). Exactly what is producing this warming warrants further investigation but it is potentially343

related to the lack of large volcanic eruptions in the projected future forcings (Fig. 1d, red).
13

[TS] anomaly (K)

2.0
1.5

(a) CESM2 (anomalies from 1920-1940)

2.5

LENS2
GHG2
AAER2
BMB2
EE2
GHG2+AAER2
GHG2+AAER2+BMB2+EE2

2.0[TS] anomaly (K)

2.5

1.0
0.5
0.0
-0.5
1920
2.5

1940

1960

1980

2000

2020

0.5

-0.5
1920

2040

(c) Greenhouse gas driven temperature anomalies

1940

1960

1980

2000

2020

2040

anomalies
0.2 (d) Anthropogenic aerosol driven temperatureLENS1XAAER1

LENS1 XGHG1
GHG2

AAER2

0

1.5
1.0
0.5

-0.2
-0.4
-0.6

0
-0.5
1920

1.0

0.0

[TS] anomaly (K)[TS] anomaly (K)

2.0

1.5

(b) CESM1 (anomalies from 1920-1940)
LENS1
LENS1 XGHG1
LENS1 XAAER1
LENS1 XBMB1
EE1 (as residual)(LENS1 XGHG1)+ (LENS1 XAAER1)

1940

1960

1980

2000

2020

-0.8
1920

2040

1940

1960

1980

2000

2020

2040

333

Fig. 2. Time evolution of global mean surface air temperature anomalies relative to the 1920 to 1940 average for (a) CESM2
and (b) CESM1 and their respective single forcing large ensembles. For CESM1, the equivalent of “Everything Else” has been
estimated as a residual (green) of the difference between LENS1 and the sum of the greenhouse gas, anthropogenic aerosol andbiomass burning aerosol contributions. The solid purple line shows the sum of the greenhouse gas and anthropogenic aerosol
contributions and in panel (a) the dashed purple line shows the sum of all 4 components (GHG2 + AAER2 + BMB2 + EE2). (c)Reproduces the greenhouse gas-forced anomalies for CESM1 and CESM2 for a more direct comparison while (d) is the same but
for the anthropogenic aerosol contribution. The shading uncertainty range is a 95% confidence interval on the ensemble meancalculated by bootstrapping members with replacement.

344

In the CESM1 single forcing large ensemble, the capacity for exploring additivity is limited,345

but we can see that the greenhouse gas contribution inferred from LENS1−XGHG1 and the346

anthropogenic aerosol contribution inferred from LENS1−XAAER1 approximately add up to the347

overall LENS1 [Ts ] anomalies (compare black and purple in Fig. 2b). This, however, is not true in348

CESM2 (compare black and solid purple in Fig. 2a) where the sum of the [Ts ] anomalies in GHG2

349

and AAER2 fall short of the LENS2 [Ts ] anomalies from the late 20th century onwards. Adding350

in the contributions from BMB2 and EE2 brings the sum a little closer to LENS2 (dashed purple

351

in Fig. 2a) but a discrepancy still exists. It is clear from comparison of the relation between the352

solid purple and black lines in Figs. 2a and b that the sum of the greenhouse gas and anthropogenic

326
327
328
329
330
331
332

14

353

aerosol contributions and how that relates to the all forcing signal differs considerably between

354

CESM1 and CESM2.

355

The prescribed GHG concentrations are rather similar between CESM1 and CESM2 over the356

period shown (supplemental Fig. 2) and a closer comparison of the GHG-forced signals between

357

CESM1 and CESM2 (Fig. 2c) reveals that the GHG-forced [Ts ] anomalies are also comparable

358

between CESM1 and CESM2. Note that this is not true locally as GHG2 warms more than359

LENS1−XGHG1 in the low latitudes but warms much less in the Northern Hemisphere (NH)360

high latitudes - a feature that will be discussed further in section 5 . CESM1 and CESM2 differ

361

considerably in the global mean [Ts ] anomalies due to anthropogenic aerosol forcing (Fig. 2d).

362

Anthropogenic aerosols continue to cool the planet out to 2050 in CESM2 while the anthropogenic

363

aerosol induced cooling in CESM1 maximizes in the 1980s and then declines, such that by 2050, the364

global mean anthropogenic aerosol-forced [Ts ] anomalies differ by about 0.4K between CESM1 and365

CESM2. This difference in anthropogenic aerosol-forced [Ts ] change could be due to differences366

in the experimental design, differences in the aerosol emissions, differences in the model physics

367

or some combination of these, as explored in the following sections.

368

4. The impact of the single forcing method on the aerosol-forced response

369

a. Global mean temperature and radiative fluxes

370

A major difference between the CESM1 and CESM2 single forcing large ensembles is the371experimental design. In CESM1, an “all-but-one” approach was used while in CESM2 an “only”

372

approach was used. To test the influence of this experimental design, we consider the additional

373

three-member ensembles: XAAER2, where the anthropogenic aerosol simulation of CESM2 was

374

performed in the same way as with CESM1; and, AAER1, where the anthropogenic aerosol

375

simulation of CESM1 was performed in the same way as CESM2.

376

Time series of 21-year running mean anomalies (chosen as a reasonable balance between reducing

377

noise, while retaining features of the time-evolution) of various global mean quantities are shown378

in Fig. 3 for each of the methods for both CESM1 and CESM2. First, it is worth noting the379

substantial differences in the global mean AOD between the CMIP5 and the CMIP6 forcings (Fig.

3803a). This is primarily due to the differences in emissions but there is also a contribution from the381

enhanced lifetime of black carbon in CESM2 (see the Appendix). [AOD] continues to rise out to
15

4
2
1940

1960

1980

2000

2020

2040

(c) TOA [SWnet ]

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

Flux (Wm 2)

Flux (Wm 2)

0

1940

1960

1980

2000

2020

0.0

Flux (Wm 2)

Flux (Wm 2)

-0.5
-1.0
-1.5
-2.0
-2.5

1940

1960

1980

2000

2020

2040

3(light) and N(dark)
member uncertainties on(LENS1 XAAER1) and AAER2
AAER1
LENS1 XAAER1AAER2
LENS2 XAAER2
1940

1980

2000

2020

2040
0.2

(e) TOA [SWnet ] [LWnet ]

0

1940
0.3
0.2
0.1
0
-0.1
-0.2
-0.3
-0.4
-0.5

1960

(d) TOA [LWnet ]

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

2040

(f) Surface [SW ](b) [TS]

Flux (Wm 2)

6

TS (K)AOD (×10 3)

8

0
-0.1
-0.2
-0.3
-0.4
-0.5
-0.6
-0.7

1960

1980

2000

2020

2040

(g) Surface [SW ]

-0.2
-0.4
-0.6

1940

1960

1980

2000

2020

2040

(h) Surface SW Albedo

0.25
0.20
Albedo (%)(a) [AOD]
10

0.15
0.01
0.05
0

1940

1960

1980

2000

2020

2040

-0.05

1940

1960

1980

2000

2020

2040

395

Fig. 3. Global mean centered 21-year running means of annual means of anomalies relative to 1920 to 1940 of various fields from
AAER1 (teal solid), LENS1−XAAER1 (teal dashed), AAER2 (maroon solid) and LENS2−XAAER2 (maroon dashed). Shaded
ranges are shown around the experiments with the large ensembles (LENS1−XAAER1 and AAER2) with the light component
showing the uncertainty on a 3-member mean of XAAER1 or AAER2 and the dark component showing the uncertainty for an
N-member mean of XAAER1 or AAER2 where N is the number of members in the anthropogenic aerosol single forcing ensemble.
(a) Aerosol optical depth at 500nm, (b) near surface air temperature, (c) top of atmosphere (TOA) net downward shortwave flux, (d)top of atmosphere net upward longwave flux, (e) top of atmosphere net downward radiative flux, (f) surface downward shortwave
flux, (g) surface upward shortwave flux, (h) albedo.

382

2050 in CESM2 but declines after about 1980 in CESM1 and the experimental design does not

383

substantially impact the [AOD] evolution (compare solid and dashed in Fig. 3a). A comparison384

of the burdens of different aerosol species in Fig. A1 and supplemental Fig. 4 indicates, perhaps

385

unsurprisingly, that the difference in [AOD] trends between AAER1 and AAER2 is dominated by386

the differing trends in anthropogenic aerosols as opposed to sea salt or dust which can respond as387

the climate changes under anthropogenic aerosol forcing.

388
389
390
391
392
393
394

396

It is clear from the [Ts ] time series in Fig. 3b that the experimental design has a substantial

397

impact on the inferred anthropogenic aerosol influence on [Ts ] in CESM2. When the aerosol

398

forcing is imposed in isolation (solid maroon in Fig. 3b), much colder temperature anomalies are
16

399

reached than when the aerosol influence is inferred from LENS2−XAAER2 (dashed maroon in

400

Fig. 3b), with the difference in the aerosol-forced cooling being ∼0.3K by 2030-2050. Note that401

the difference between LENS2−XAAER2 and AAER2 is much greater than would be expected402

due to sampling uncertainty alone (compare dashed maroon with the light uncertainty shading on403

the solid maroon in Fig. 3b). The result is that when the anthropogenic aerosol influence is inferred404

from LENS2−XAAER2, the overall cooling is more comparable to the cooling found in CESM1

405(compare maroon dashed with teal in Fig. 3b). This indicates a strong sensitivity of the inferred406

aerosol cooling to the experimental design (“all-but-one” versus “only”) in CESM2. In CESM1,

407

the AAER1 simulation (solid teal) is also significantly cooler than LENS1−XAAER1 (dashed teal)

408

toward the end of the simulation, but the impact of the experimental design in CESM1 is relatively

409

minor compared to that found in CESM2. We will revisit this difference between CESM1 and410

CESM2 in section 5 and for now focus on the dependence on the experimental design within

411

CESM2.

412

In Fig. 4a it can be assessed how [Ts ] varies as a function of [AOD]. Here, we only show the413

variations over the time period when [AOD] is increasing, which means for CESM1, we are showing

414

out to the 21-year mean centered on 1984, while for CESM2 we are showing out to the end of the415

simulation. Fig. 4a shows that up to [AOD] anomalies of about 8×10−3 (the maximum in CESM1)

416

the evolution of [Ts ] as a function of [AOD] is rather similar in each of AAER1, LENS1−XAAER1

417

and LENS2−XAAER2. In contrast, AAER2 cools a lot more at a given [AOD] for [AOD] greater

418

than around 4×10−3 . In addition, the evolution of [Ts ] as a function of [AOD] is non-linear,419

particularly in LENS2−XAAER2. In the 21st century, the cooling in LENS2−XAAER2 levels off

420

and then turns around and the planet starts to warm even while the [AOD] continues to increase

421(Fig. 3a versus b and Fig. 4a).

427

To begin to understand the difference between AAER2 and LENS2−XAAER2, consider the428

top of atmosphere (TOA, although actually here the fluxes used are at the model top) radiative429

fluxes and their imbalance shown in Figs. 3c-e. Fig. 3e shows the difference between the net TOA

430

downward shortwave radiation ([SWnet ↓]) and the net TOA upward longwave radiation ([LWnet ↑]),

431

i.e., the TOA radiative imbalance. Throughout we refer to fields that are positive when downward

432

with the down arrow (↓) and fields that are positive when upward with the up arrow (↑). The way

433

in which the overall TOA radiative imbalance evolves is similar between the methods (compare

17

AAER2 (<1984)
AAER2 ( 1984)
LENS2 XAAER2 (<1984)
LENS2 XAAER2 ( 1984)
AAER1
LENS1 XAAER1

0 1 2 3 4 5 6 7 8 9 10
[AOD] (×10 3)(d) TOA [SWnet ] vs [TS]
[LWnet ] (Wm 2)

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

(b) TOA [SWnet ] vs [AOD]-0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0[TS] (K)

0.25

0.20

0.20

0.15
0.10
0.05
0

423
424
425
426

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

(c) TOA [LWnet ] vs [AOD]

0 1 2 3 4 5 6 7 8 9 10
[AOD] (×10 3)(e) TOA [LWnet ] vs [TS] (K)-0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0[TS] (K)(g) Surface SW Albedo vs TS

0.15
0.10
0.05
0

0 1 2 3 4 5 6 7 8 9 10
[AOD] (×10 3)
422

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

0.25 (f) Surface SW Albedo vs AOD
Surface SW Albedo (%)

Surface SW Albedo (%)[SWnet ] (Wm 2)

0 1 2 3 4 5 6 7 8 9 10
[AOD] (×10 3)

0
-0.2
-0.4
-0.6
-0.8
-1.0
-1.2
-1.4
-1.6

[LWnet] (Wm 2)(a) [TS] vs [AOD][SWnet] (Wm 2)[TS] (K)

0
-0.1
-0.2
-0.3
-0.4
-0.5
-0.6

-0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0
TS (K)

Fig. 4. Global mean 21-year running means of annual means. Both CESM1 and CESM2 are only shown for the period overwhich [AOD] is continuing to increase, which means the full record is shown for CESM2, but only up to the 21-year mean centeredon 1984 is shown for CESM1. The CESM2 points transition from circles to triangles after 1984. (a)-(c) show [Ts ], TOA net
downward shortwave and TOA net upward longwave, respectively, versus [AOD]. (d) and (e) show TOA net downward shortwave
and TOA net upward longwave versus [Ts ] and (f) and (g) show surface shortwave albedo versus (f) [AOD] and (g) [Ts ].

434

solid and dashed in Fig. 3e). However, the evolution of the separate components [SWnet ↓] (Fig.

4353c) and [LWnet ↑] (Fig. 3d) is not - they reveal that this same TOA radiative imbalance is achieved

436

for rather different reasons.
18

437

Until about 1980, aerosol forcing causes TOA [SWnet ↓] to decline as aerosols and associated438

cloud changes reflect more shortwave radiation back to space and as surface albedo increases. In

439

AAER2, after about 1970, there is a greater decline in TOA [SWnet ↓] than in LENS2−XAAER2,

440

but this difference is opposed by a greater decline in [LWnet ↑], such that the overall TOA radiative441

imbalance ends up roughly the same in each case. The greater decline in [LWnet ↑] in AAER2

442

can be attributed to the greater decline in [Ts ]; a colder planet emits less longwave radiation to443

space. Indeed, Fig. 4e shows that [LWnet ↑] depends on [Ts ] in a rather similar way in AAER2 and444

in LENS2−XAAER2 but AAER2 cools further and, in association with this, [LWnet ↑] declines445

more. So, overall, while the TOA radiative imbalance evolves in a similar way in AAER2

446

and in LENS2−XAAER2, it reflects different quasi-equilibria with differing balances between

447

TOA [LWnet ↑] and TOA [SWnet ↓]1. Overall, we infer that the reason AAER2 cools more than

448

LENS2−XAAER2 lies in the behavior of the shortwave radiation. TOA [SWnet ↓] declines more449

in AAER2 and this is balanced by a greater cooling and greater reduction in TOA [LWnet ↑].

450

The greater decline in TOA [SWnet ↓] in AAER2 compared to LENS2−XAAER2 must arise

451

from a difference in the shortwave radiation being reflected back to space, either from within452

the atmosphere or from the surface. Fig. 3 demonstrates that it is the difference reflected back453

from the surface that is key to the TOA [SWnet ↓] differences.

454

the extent to which shortwave radiation is reflected back to space from within the atmosphere,

455

either by aerosols themselves or the associated cloud radiative effects, then we would expect to see

456

differences in the surface downward shortwave radiation ([SW ↓]), but Fig. 3f makes it clear that457

the anomalies in [SW ↓] are very similar between AAER2 and LENS2−XAAER2. The decline458

in surface [SW ↓] is independent of which method is used in both CESM1 and CESM2, although459

the response differs substantially between the models in association with their differing aerosol

460

forcings. In contrast, Fig. 3g demonstrates a clear difference in the surface upward shortwave

461([SW ↑]) between AAER2 and LENS2−XAAER2 in association with a difference in their surface

462

shortwave albedo responses (Fig. 3h). In AAER2, the surface shortwave albedo increases much463

more than in LENS−XAAER2. As a result, surface [SW ↑] stays roughly constant in AAER2

If the key were differences in

1An aside is that in Fig. 3d, LENS2−XAAER2 is the odd one out, with AAER1, LENS1−XAAER1 and AAER2 all exhibiting similar changes
in [LWnet ↑]. This may appear at odds with the fact that it is AAER2 that exhibits a different temperature response (Fig. 3b). The reason why
AAER1 and LENS1−XAAER1 exhibit a greater decline in TOA LWnet ↑ than LENS2−XAAER2, even though their temperature responses are
similar is actually because of cloud longwave radiative effects. Examination of clear sky LWnet ↑ (supplemental Fig. 5) reveals what we expect:
AAER1, LENS1−AAER1 and LENS2−XAAER2 which all cool less than AAER2, also exhibit a smaller decline in clear sky [LWnet ↑].

19

464

(Fig. 3g), even as surface [SW ↓] declines, because a larger proportion is being reflected back to465

the atmosphere, and ultimately to space.

466

The difference in surface [SW ↑] between AAER2 and LENS2−XAAER2 (solid in Fig. 5a)467

explains most of the difference in TOA clear sky (dashed in Fig. 5a) and TOA all sky [SWnet ↑]468(dotted in Fig. 5a). Consideration of how the surface shortwave albedo varies as a function of469[AOD] (Fig. 4f) and [Ts ] (Fig. 4g) reveals that there is a systematic difference between the “all-470

but-one” and “only” approach in both CESM1 and CESM2. For [Ts ] anomalies down to around471-0.3K and [AOD] anomalies up to about 5×10−3 the albedo increases more for the “only” approach472

than for the “all-but-one” approach. But then, beyond that, the difference in behavior of the albedo

473

between AAER2 and LENS2−XAAER2 increases rather dramatically. In the late 20th century,474

the surface shortwave albedo continues to increase in AAER2 (Fig. 3h), but in LENS2−XAAER2

475

the albedo increase levels off at a much lower value and then starts to decline. This is apparent as476

a rather dramatic difference between AAER2 and LENS2−XAAER2 in the relationship between

477

surface shortwave albedo and both [AOD] and [Ts ] (Figs. 4 f and g). In AAER2, as the planet478

cools, surface shortwave albedo keeps on increasing, presumably providing a positive feedback479

onto the cooling (Fig. 4g). In LENS2−XAAER2, as the planet cools, albedo also increases, but

480

to a lesser extent, and then in the 1990s, albedo starts to decline and the planet begins to warm up

481

again, even though the [AOD] has continued to increase.

488

Increased surface shortwave albedo generates cooler temperatures and vice versa, so separating489

out cause and effect is challenging in these quasi-equilibrium experiments where the system has

490

adjusted to a new balance. Nevertheless, given that there is no evidence that the origins of the491

different [Ts ] response between AAER2 and LENS2−XAAER2 lies in differences in how the492

incoming surface shortwave radiation behaves, we posit that it lies in non-linearity in surface

493

shortwave albedo feedbacks and that there are two components that contribute to this shortwave

494

albedo non-linearity, which will now be discussed: (1) a non-linearity in both snow and sea ice

495

albedo feedbacks and; (2) a non-linearity related to the North Atlantic ocean circulation which496

leads to differing northward heat transport into the Arctic and associated differences in high497

latitude albedo. In the following sub-sections we provide the evidence for both of these sources of498

non-linearity in surface shortwave albedo.

20

Spatially inetegrated SW (×1014W)

2.5
2.0
1.5

(a) AAER2 (LENS2 XAAER2), global and contributions
TOA [SWnet ]TOA clear sky [SWnet ]Surface [SW ]50 S-50 N
50 S-90 S
50 N-90 N

1
0.5
0
1940

Spatially integrated SW (×1014W)

1.4
1.2

1980

2000

2020

2040

(b) AAER2 (LENS2 XAAER2), 50 N-90 N
Surface [SW ] 50 N-90 N
Land
Not Land

1.0
0.8
0.6
0.4
0.2
0
1940

Spatially integrated SW (×1013W)

1960

4.0
3.5
3.0
2.5
2.0
1.5
1.0
0.5
0
-0.5

1960

1980

2000

2020

2040

(c) AAER2 (LENS2 XAAER2), 50 S-90 S
Surface [SW ] 50 S-90 S
Land
Not Land

1940

1960

1980

2000

2020

2040

487

Fig. 5. Difference between AAER2 and LENS2−XAAER2 in the spatially integrated shortwave fluxes. (a) shows the global
integral with black dotted showing the top of atmosphere net upward shortwave radiation and black dashed showing its clear sky
component. Black solid shows the surface upward shortwave radiation and it is further divided into contributions from different
latitude bands. (b) Surface upward shortwave radiation spatially integrated from 50◦ N to 90◦ N and the contributions from land
regions and regions that are not land. (c) is as (b) but for 50◦ S to 90◦ S. Note the differing magnitudes covered by the y axes in (b)and (c).

499

b. Albedo non-linearities

5001) Snow and sea ice non-linearities

482
483
484
485
486

501

Figure 5a shows the spatially integrated difference in the surface upward shortwave radiation,502

21 (solid black) and this difference is further de[SW ↑] between AAER2 and LENS2−XAAER2

503

composed into the contribution from different latitude bands. The high latitudes (poleward of 50◦

504

latitude) clearly dominate in this difference. Both SH (salmon) and NH (blue) play a role, but the505

NH difference in surface SW ↑ dominates in the late 20th and 21st centuries. This difference in506

surface SW ↑ between AAER2 and LENS2−XAAER2 in the high latitudes of each hemisphere can507

be further decomposed into the contribution from land regions and the contribution from “not land”

508

regions, i.e., regions that are either ocean or sea ice. In the NH, the land and not land contribute509

roughly equally to the difference between AAER2 and LENS2−XAAER2 from the 2000s onward

510(Fig. 5b) and the difference over land regions dominates prior to that. In the SH, the difference

511

over ocean and sea ice regions dominates (Fig. 5c).

512

A seasonal decomposition of the difference in surface SW ↑ integrated over 50◦ N to 90◦ N between

513

AAER2 and LENS2−XAAER2 indicates the summer season as dominating in the SW ↑ difference514

initially (supplemental Fig. 6), which makes sense given that this is when there is the greatest

515

incident shortwave radiation which can then be affected by differences in surface shortwave albedo

516

feedbacks. Later in the simulation, differences in SW ↑ between the methods becomes increasingly

517

important in the shoulder seasons as well. In order to understand the origins of this difference in518

surface shortwave albedo behavior we now focus on the NH during summer (June-July-August,

519

JJA).

520

Consider the time series of 21-year running mean JJA 50◦ N-90◦ N Ts shown in Fig. 6a. This

521

shows that NH high latitude temperature declines at a similar rate in AAER2 and LENS2−XAAER2

522

until around 1960-1980, at which point LENS2−XAAER2 starts to warm, while AAER2 does not.

523

We consider the behavior of surface SW↑, snow cover and sea ice cover during the 1960-1980

524

average in an attempt to examine their differences before subsequent feedbacks associated with525

the differing Ts response are present. For 50◦ N-90◦ N average surface SW ↑, a difference between

526

AAER2 and LENS2−XAAER2 is already apparent during 1960-1980, even though a difference in527

Ts is not (see Fig. 6b-d for local temperature changes). It is clear from Figs. 6 f-h that in AAER2,

528

there is a greater enhancement in surface SW ↑ around the sea ice edge and over high latitude land529

regions than in LENS2−XAAER2. There is also less of an increase in surface SW ↑ in the interior

530

regions of the sea ice in AAER2 compared to LENS2−XAAER2, leading to a difference in SW ↑531

between the methods that is negative there. These differences in SW↑ correspond reasonably well532

to differences in sea ice cover between AAER2 and LENS2−XAAER2 (Fig. 6n-p).

22

(d) (b) - (c)

2

-0.2

1

-0.4

0

-0.6

1

-0.8

2

1940 1960 1980 2000 2020 2040

(g) LENS2 XAAER2(f) AAER2(h) (f) - (g)10SW (Wm 2)

3
2

0

1
0
-1

Snow cover (%)

41940 1960 1980 2000 2020 2040

(i) Snow cover, 50 N-90 N(k) LENS2 XAAER2(j) AAER2

10

2

0

1

5

10

1940 1960 1980 2000 2020 2040(m) Ice cover (%), 50 N-90 N(o) LENS2 XAAER2(n) AAER2(p) (n) - (o)

10

3

5

2

0

1

5

0

10

1940 1960 1980 2000 2020 2040
533(l) (j) - (k)54
Ice cover (%)

10

3

0

SW (Wm 2)(e) Surface SW , 50 N-90 N
4Snow cover (%)-1.0

Ice cover (%)TS(c) LENS2 XAAER2(b) AAER2

AAER2
LENS2 XAAER2TS (K)(a) TS, 50 N-90 N
0Fig. 6. Aerosol influence on the JJA season from 50◦ N to 90◦ N. Left column shows time series of 21-year running means543

for (solid) the ensemble mean of AAER2 and (dashed) the ensemble mean of LENS2-XAAER2 anomalies from the 1920-1940
average. The dark and light shadings around the AAER2 line show the 95% confidence interval using (light shading) 3 membersand (dark shading) 15 members. The dashed vertical line depicts the 21-year mean centered on 1970, i.e., the 1960-1980 averageshown in the right three columns. (2nd column) shows 1960-1980 anomalies compared to 1920-1940 for AAER2, (3rd column)
is as (2nd column) but for LENS2−XAAER2 and (4th column) shows the difference in the anomalies between AAER2 andLENS2−XAAER2. (top) Ts , (2nd row) surface upward shortwave, (3rd row) grid cell area covered by snow in percent, (4th row)
grid cell area covered in sea ice in percent. In (i) and (m) the 50◦ N-90◦ N average is taken only over land grid points and grid points
that are not land, respectively. Blue points in (l) show the grid points used for the analysis in Fig. 7 ( (286◦ E,58◦ N), (173◦ E,68◦ N),(95◦ E,73◦ N) ). The blue and pink contours in (n) and (o) show the 80% sea ice contours for 1920-1940 of AAER2 and 1960-1980
of XAAER2, respectively.

544

Over the NH land regions surrounding the Arctic, there is a greater increase in summer snow

545

cover in AAER2 than in LENS2−XAAER2 (Fig. 6l) corresponding well to regions where the546

difference in surface SW ↑ is largest over land. Figure 7 demonstrates that snow cover fraction

534
535
536
537
538
539
540
541
542

23

547

depends non-linearly on local temperature at three representative locations, given by the blue548

circles in Fig. 6l: a grid point to the east of Hudson Bay; one in Eastern Siberia; and one to the549

south of the Kara Sea. Specifically, a cooling/warming that occurs at a lower temperature would

550

be associated with a larger increase/decrease in snow cover than if that same cooling/warming

551

occured at a higher temperature. This can be understood as a result of snow cover being bounded552

by zero. For warmer temperatures, there is a higher probability of there not being snow on the553

ground (green stars in Fig. 7) and for times when there is no snow on the ground a further warming554

can no longer lead to a change in snow cover, leading to a weaker dependence of seasonal average

555

snow cover on temperature at warmer temperatures. The arrows in Fig. 7 help to illustrate this556

non-linearity by showing the 1920-1940 average at the start point of the arrow and the 1960-1980

557

average at the end point of the arrow. For the change that is inferred to be due to the aerosol forcing,558

the cooling in AAER2 which has a colder starting point, leads to a proportionately larger increase559

in snow cover than the cooling in LENS2 does and both lead to a proportionately larger change in

560

snow cover than the warming in XAAER2, which is warming rather than cooling. The result is561that the magnitude of the increase in snow cover in AAER2 (given by the length of the blue arrow)562

is proportionately larger for the temperature change than that in LENS2−XAAER2 (given by the563

sum of the black and the pink arrow lengths). This effect is likely what dominates prior to 1980

564

in the differences seen between the methods in the NH, given the dominance of land regions in

565

contributing to the SW↑ differences (Fig. 5b).

576

For sea ice, both AAER2 and LENS2−XAAER2 exhibit an increase in sea ice cover, but they577

do so at different locations - the increase in AAER2 (Fig. 6n) is generally at lower latitudes than578

the increase in LENS2−XAAER2 (Fig. 6o). The reason for this is fairly straightforward - in

579

AAER2, the cooling is occurring relative to a cold climate (a pre-industrial climate which has then580

cooled slighly under aerosol forcing out to the baseline 1920 to 1940 period), while in LENS2

581

compared to XAAER2, the aerosol influence is felt relative to a climate in which sea ice has been582

influenced by greenhouse gas driven warming. To illustrate the differing sea ice fractions between

583

the different baseline climates that the aerosol influence is being compared against in Fig. 6, we584

show the 80% sea ice fraction contours for the 1920-1940 climate of AAER2 in blue and for the5851960 to 1980 climate of XAAER2 in pink in panels n and o. The 1920-1940 climate is the baseline

586

for AAER2. For the “all-but-one” method, the baseline climate is more complicated but since

24

20

0
272 274 276 278 280 282 284 286 288 290
Local TS

LENS2
AAER2
XAAER2

40
20

60
40
20

0

0
272 274 276 278 280 282 284 286 288 290
Local TS

80

60
Snow cover (%)

40

80

60

100

80

40
20
0

LENS2
AAER2
XAAER260
40
20

Probability of snow free days (%)060(c) 95 E, 73 N

100

80

Snow cover (%)

Snow cover (%)

20

LENS2
AAER2
XAAER2

Probability of snow free days (%)

80

60
40

(b) 173 E, 68 N

100

Probability of snow free days (%)(a) 286 E, 58 N
80

0
272 274 276 278 280 282 284 286 288 290
Local TS

575

Fig. 7. Evolution of snow cover as a function of local Ts for the three points shown in Fig. 6l: (a) 286◦ E, 58◦ N (east of Hudson
Bay); (b) 173◦ E, 68◦ N (eastern Siberia); and, (c) 95◦ E, 73◦ N (south of the Kara Sea). Dots show the JJA seasonal mean percentageof the grid point covered by snow versus local Ts . Small dots show the individual seasons for all members and all years and large
dots show the ensemble means for each year. Blue = AAER2, pink = XAAER2 and black = LENS2 (LENS2 is shown out to 2100).
The green stars (right axis) show the probability of snow free days in the JJA season assessed for each year by pooling together all
members from either LENS2 or AAER2 (XAAER2 is not shown for this metric given its smaller ensemble size). The start pointof each of the arrows shows the ensemble mean 1920-1940 average value and the end point shows the ensemble mean 1960-1980
average. The overall change in AAER2 is simply quantified by the length and direction of the blue arrow while the magnitude of thechange in LENS2−XAAER2 is given by summing up the length of the black and pink arrows when they are in opposite directions,as in all cases here and the direction of change is that of the black arrow.

587

we are comparing the 1960-1980 minus 1920-1940 anomalies of LENS2 with that in XAAER2,

588

and the sea ice fraction doesn’t differ substantially between LENS2 and XAAER2 in 1920-1940

589(not shown), the baseline for the aerosol influence in LENS2−XAAER2 is effectively 1960-1980

590

of XAAER2. The 80% sea ice fraction metric indicates that a high sea ice fraction is present

591

for a wider latitude range in the colder AAER2 climate and, as a result, the additional growth592

due to the anthropogenic aerosol influence occurs at lower latitudes compared to the growth that593

occurs in LENS2 relative to XAAER2. Note that the difference in ice fractions between 1920-1940

594

of AAER2 and 1960-1980 of XAAER2 is also accompanied by differences in sea ice thickness

595(thicker ice in the central Arctic in AAER2) which will also impact on where additional sea ice

596

under aerosol-forced cooling will grow. The overall result is that AAER2 gains more sea ice at597

low latitudes and less sea ice at high latitudes compared to LENS2−XAAER2 (Fig. 6p). While598

the latitude at which sea ice is gained under aerosol forcing clearly represents a state dependence599

leading to differences between the methods, its effects on SW↑ during the summer months is likely

600

small, given that during the summer, the latitudinal gradients in incoming shortwave are small.

601

But this could be a contributor to the method dependence of SW↑ anomalies during the spring602

and autumn, when insolation is relatively greater at lower latitudes. In addition to this difference603

in latitude at which sea ice is gained, supplemental Fig. 7 shows the dependence of JJA averaged

566
567
568
569
570
571
572
573
574

25

604

sea ice fraction on temperature, in a similar manner to Fig. 7 for snow, where it can be seen that a

605

non-linearity is also present in association with sea ice being bounded by zero and, therefore, the606

probability of having days with zero sea ice increases with increasing temperature.

607

In the SH, summer sea ice fraction differences and their relation to differences in surface upward608

shortwave radiation are also clear. AAER2 shows greater increases in sea ice cover overall with a609

strong correspondence between regions where sea ice has increased more and regions where the610

increase in upward shortwave from the surface is greater (supplemental Fig. 8).

611

In summary, the behavior of snow cover and sea ice and their influence on surface upward612

shortwave radiation appear to be state dependent, i.e., it matters whether aerosol forcing is imposed

613

within a cold pre-industrial climate or whether it is imposed within a climate state that has also614

experienced greenhouse gas forcing. An aerosol cooling that occurs at a colder temperature, as in

615

AAER2, increases the average snow cover more than an aerosol induced cooling that occurs within616

a planet that has warmed under greenhouse gas forcing, as in LENS2 compared to XAAER2.

617

The same is true for sea ice and, in addition, sea ice that is gained during the cooling of a colder

618

climate tends to occur at lower latitudes than the sea ice gains that occur during the cooling of a619

warmer climate. These snow cover and sea ice non-linearities lead to an overall larger influence620

on shortwave radiation for the colder base state in AAER2.

6212) The North Atlantic Ocean circulation

622

The annual mean Ts response to anthropogenic aerosol forcing by 2030-2050 can be seen in Fig.

623

8. By this time period, AAER2 (Fig. 8a) is colder than LENS2−XAAER2 (Fig. 8b) over much624

of the globe, with the largest differences found in the NH high latitudes, over continental regions,625

around the margins of Antarctic sea ice and in the tropical and sub-tropical Pacific (Fig. 8c). The

626

fact that the Ts differences extend beyond the high latitudes is not inconsistent with the important627

role for high latitude feedbacks in producing them because similar anomalies, with opposite sign628

relative to those in Fig. 8c in the tropical Pacific have been found in response to sea ice loss (as629

opposed to gain in our case) (Deser et al. 2015) with an important role for ocean dynamics in

630

transferring the signal there (Wang et al. 2018).
26

(a) AAER2(b) LENS2 XAAER2(c) (a) (b)(d) AAER1(e) LENS1 XAAER1(f) (d) (e)

2.0

1.5

1.0

0.5

0.0

TS anomaly (K)

0.5

1.0

1.5

2.0

640

Fig. 8. Differences in annual mean Ts between 2030-2050 and 1920-1940. (a)-(c) CESM2 AAER2, LENS2−XAAER2 andthe difference between them. (d)-(f) as (a)-(c) but for CESM1. Stippling indicates anomalies that are not statistically significant at
the 95% confidence level.

631

Consideration of the differences between AAER2 and LENS2−XAAER2 in 2030-2050, however,632

reveals another important feature: there is a clear difference in the sub-polar North Atlantic with a633

substantial warm anomaly in the region south of Greenland in LENS2−XAAER2 (Fig. 8b), which634

is much less apparent in AAER2 (Fig. 8a). In LENS2−XAAER2, there is also a warm anomaly

635

over much of the Arctic. It is clear that there are large differences between the two methods in the636

Ts response over the NH high latitudes and in the sub-polar gyre region to the south of Greenland,637

in particular (Fig. 8c), while the same method-dependence is not found in CESM1 (Fig. 8d-f).

638
639

641

The warm sub-polar gyre Ts anomaly in LENS2−XAAER2 to the south of Greenland resembles

642

what would be expected from a strengthening of the AMOC (e.g., Delworth et al. 2017). Indeed,

643

consideration of the AMOC response, defined as the change in the magnitude of the maximum

644

meridional overturning streamfunction at 45◦ N below 500 m depth, reveals a strong dependency

645

of the aerosol-forced AMOC changes on the experimental design in CESM2 (Fig. 9c).

646First, it is worth considering how the NH aerosol forcing evolves as this is likely to be more647

directly connected to forcing of AMOC changes than the global mean aerosol evolution. Figure648

9a shows that the 50◦ N to 90◦ N AOD in CESM2 increases to a maximum in the 1970s and 1980s

649

and then declines but levels off at higher AOD values than the 1920-1940 period for the remainder
27

650

of the simulation. In AAER2, the AMOC increases in strength to a maximum in the 1970s and651then declines (Fig. 9c, solid), somewhat following the NH high latitude aerosol forcing, except it

652

returns to the baseline AMOC strength despite the positive AOD anomalies in the 21st century.

653

In contrast, the AMOC in LENS2−XAAER2 increases in strength more rapidly and then, rather

654

than decreasing with the decline in aerosol forcing, it plateaus and even slightly increases out to655

the end of the simulation (Fig. 9c, dashed). Consideration of XAAER2 in isolation (pink in Fig.

6569e) reveals that as the planet warms in the absence of aerosol forcing and the presence of rising657

greenhouse gases, the AMOC in CESM2 starts to decline rapidly around 1980. This is rather

658

similar to what is seen in GHG2 (red in Fig. 9g). In contrast, in LENS2 when all forcings are659

present, the anthropogenic aerosol forcing seems to dominate and acts to strengthen the AMOC

660

until about 1980, delaying this rapid greenhouse gas-forced decline in AMOC until later in the661

simulation (black in Figs. 9e and g). The result is a non-linear behavior of AMOC in the CESM2

662

simulations with the sum of the individual forcing contributions not adding up to the LENS2

663

response (compare black and dashed purple in Fig. 9g).

675

The LENS2 AMOC decline begins around 1980 but it does not fall below the 1920-1940 baseline

676

until about 2000 and at that point it declines at a rather similar rate to what was seen earlier in

677

XAAER2. The AMOC in XAAER2 starts to decline very rapidly around 1980, and because the678

aerosols in LENS2 delay the onset of this rapid decline compared to XAAER2, the aerosol influence679

inferred from LENS2−XAAER2 is an apparently greater strengthening of the AMOC than that680

inferred when the aerosols are imposed in isolation in AAER2. It is not that the aerosol forcing

681

by itself produces the strengthening inferred from LENS2−XAAER2, but it is that it staves off the682

rapid greenhouse gas-forced AMOC decline. This can explain the warm anomalies in the NH high683

latitudes and the sub-polar North Atlantic due to aerosol forcing estimated from LENS2−XAAER

684

in Fig. 8b, whereas in AAER2, the AMOC strengthening is weaker and declines after the 1980s.

685

The increased AMOC strength in LENS2−XAAER2 leads to enhanced northward ocean heat

686

transport into the NH high latitudes (not shown), which is a further boost to the disparity in surface

687

albedo between the two methods through effects on sea ice and snow melt.

688

It is challenging to truly isolate the relative importance of AMOC versus the other albedo non-689

linearities described above to the method-dependence of the global surface upward shortwave

690

radiation, but we can at least obtain a rough estimate of the order of magnitude of AMOC’s

28

15

(a) AOD anomalies (50 N-90 N), CESM2

1510AOD (×10 3)

AOD (×10 3)

10
5
0
-5

AAER2
LENS2 XAAER2AMOC (Sv)

2

1960

1980

2000

2020

AMOC (Sv)

AAER2
LENS2 XAAER2

1960

1980

2000

2020

(e) LENS2, AAER2 and XAAER2 AMOC

2040

3.5
3.0
2.5
2.0
1.5
1.0
0.5
0.0
-0.5
-1.0
2
0

-2

-2

-4
-6
-10
2

LENS2
AAER2
XAAER2

1940

1980

2000

2020

(g) LENS2 AMOC change contributions

2040

-8
-10

665
666
667
668
669
670
671
672
673
674

AMOC (Sv)

AMOC (Sv)

-6

2000

2020

2040

1960

1980

2000

2020

2040

1960

1980

2000

2020

2040

2020

2040

AAER1
LENS1 XAAER1

1940(f) LENS1, AAER1 and XAAER1 AMOC

-10

LENS1
AAER1
XAAER11940(h) LENS1 AMOC change contributions

0

-2
-4

1980

-6

2

0

1960

(d) CESM1 AMOC due to aerosols (45 N)

-4
-8

1960

AAER1
LENS1 XAAER1

1940

0

-8

664

0

2040

(c) CESM2 AMOC due to aerosols (45 N)

1940

5

-5

AMOC (Sv)

AMOC (Sv)

1940

3.5
3.0
2.5
2.0
1.5
1.0
0.5
0.0
-0.5
-1.0

(b) AOD anomalies (50 N-90 N), CESM1

LENS2
AAER2
GHG2
BMB2
EE2
GHG2+AAER2
GHG2+AAER2+BMB2+EE2

1940

1960

1980

-2
-4
-6
-8

2000

2020

2040

-10

LENS1
LENS1 XAAER1
LENS1 XGHG1
LENS1 XBMB1(LENS1 GHG1) + (LENS1 AAER1)

1940

1960

1980

2000

Fig. 9. (a) and (b) 21-year running mean AOD anomalies (relative to 1920 to 1940) averaged from 50◦ N to 90◦ N for CESM2 andCESM1 anthropogenic aerosol forcing, respectively. The remaining panels show 21-year running mean AMOC anomalies (relative
to 1920 to 1940) where AMOC is defined as the magnitude of the maximum meridional overturning streamfunction below 500m
depth in the North Atlantic at 45◦ N. (c) and (d) show the inferred AMOC changes due to anthropogenic aerosol forcing using bothmethods for CESM2 and CESM1, respectively, with 95% confidence intervals provided for AAER2 in (c) and LENS1−XAAER2in (d) using 3 members (light) and the number of members in the large ensemble (dark). (e) and (f) show the AMOC anomalies
for all forcings, the AAER simulation and the XAAER simulation for CESM2 and CESM1, respectively. (g) and (h) show thedecomposition of the overall change in AMOC in the large ensemble into the contributions that are inferred to be due to individual
forcings for CESM2 and CESM1, respectively. 95% confidence intervals on the ensemble means are shown in (e)-(h). For the 3
member XAAER2 and AAER1 ensembles in (e) and (f), respectively, the 95% confidence interval is calculated by bootstrapping
the pre-industrial control.

29

691

impacts by considering the association of globally integrated surface upward shortwave radiation

692

with AMOC variability in the CESM2 pre-industrial control simulation. Figure 10a shows the693

regression of 21-year running mean globally integrated surface SW↑ anomalies onto AMOC (after

694

linearly detrending to remove the pre-industrial control drift). This shows that following an increase695

in AMOC of 1Sv, the globally integrated surface SW↑ declines by just under 2.5×1013 W about696

5 years later. We then use this relationship between AMOC and globally integrated surface SW↑

697

with a lag of 5 years to construct the influence of the AMOC anomalies in each experiment on

698

globally integrated surface SW↑ and the difference in this between AAER2 and LENS2−XAAER2.

699

A comparison of this constructed AMOC influence on surface SW↑ with the difference in globally700

integrated surface SW↑ between AAER2 and LENS2−XAAER2 suggests that the influence of701

the method-dependence on AMOC can explain a little under half of the influence of the method-702

dependence on the globally integrated surface upward SW↑, with presumably the other albedo

703

effects described above contributing to the remainder (Fig. 10b). This assumes that we can linearly704

relate globally integrated SW↑ to AMOC variability and that there is no dependence of AMOC-

705

related surface upward shortwave variability on the climate base state. An analysis of similar706

regressions to those in Fig. 10a throughout the transient LENS2 simulations suggests that it is,707

indeed, a reasonable approximation to assume that the piControl regression of SW↑ onto AMOC

708

is representative of that over the 20th and early 21st centuries (Supplemental Fig. 9).

722

A variety of processes can force an AMOC decline under climate change, including reduced

723

sensible heat loss from the ocean in the presence of a warmer atmosphere (Weaver et al. 2007;

724

Brodeau and Koenigk 2016), altered freshwater forcing as precipitation and evaporation patterns

725change (Manabe and Stouffer 1993; Dixon et al. 1999) and altered lateral transports of freshwater726

into regions of deep convection as a result of sea ice loss (or melting of the Greenland Ice Sheet,

727

although Greenland ice sheet melt is not represented in these CESM2 simulations) (Jahn and

728

Holland 2013; Yang et al. 2016; Li et al. 2021). Once an AMOC decline has been induced, positive729

feedbacks, particularly from the reduced advection of salty water from southern latitudes, can730

further enhance the AMOC decline. Such feedbacks have also been argued recently by Hassan

731

et al. (2021) and Robson et al. (2022) to be important in the aerosol-forced strengthening of AMOC.

732

A more detailed analysis of the reasons behind the substantial AMOC decline under GHG forcing

733

in CESM2 is warranted and, while we leave this for future work, we provide a cursory assessment

30

Surface SW (×1013W/Sv)(a) Regression of SW onto AMOC

1.0
0.5
0.0
-0.5
-1.0
-1.5
-2.0
-2.5 SW leads

Minimum value

SW lags

Surface SW (×1014W)-40 -30 -20 -10 0 10
Lag (years)
2.00
1.75
1.50
1.25
1.00
0.75
0.50
0.25
0
-0.25

20

30

40

(b) Globally integrated surface SW
AAER2 (LENS2 XAAER2)
Constructed AMOC influence

1940

1960

1980

2000

2020

2040

721

Fig. 10. (a) Regression of 21-year running mean globally integrated surface upward shortwave radiation onto 21-year runningmean AMOC strength, defined as the maximum streamfunction below 500 m at 45◦ N, in the CESM2 pre-industrial control. The
gray shaded range shows a 95% confidence interval determined by bootstrapping with replacement 200-year segments of the
pre-industrial control, concatenating them to obtain 1000 time series of equivalent length to the 1600-year pre-industrial control
simulation, recalculating the regression and obtaining the 2.5th to 97.5th percentile range. The red point marks the minimum
value of this regression curve. (b) Black shows the time series of the difference between AAER2 and LENS−XAAER2 21-year
running mean globally integrated surface upward shortwave (reproduced from Fig. 5a) and red shows the estimated influence of thedifference in AMOC changes between AAER2 and LENS2−XAAER2 by constructing the AMOC influence on globally integrated
surface SW↑ by -2.26e13*AMOC(t-5) where -2.26e13 is the minimum regression coefficient in panel (a) and occurs at a lag of 5
years and AMOC refers to the annual mean AMOC strength anomalies relative to 1920-1940 using the maximum streamfunction
below 500m at 45◦ N definition for AMOC strength. The uncertainty range on this construction is determined by recalculating theconstruction using the bootstrapped minimum regression coefficients and lags that were used to determine the confidence interval
in panel (a) and obtaining the 2.5th to 97.5th percentile range of these bootstrapped constructions.

734

of the different forcing factors that could lead to an AMOC decline for the XAAER2 simulation735

in supplemental Fig. 10 to shed some light on the possible causes of the AMOC declines shown736

in Fig. 9e and g. This suggests that the freshwater input associated with sea ice loss is the most737

likely candidate forcing of the AMOC decline. The annual mean sea ice thickness anomalies are738

also shown in Fig. 11b and this shows that the Arctic sea ice thickness declines occur earlier in739

the absence of aerosol forcing (compare pink and black lines in Fig. 11b). This freshwater forcing

740

likely leads to a decline in the near surface density of sea water (ρ) in the Labrador Sea through

709
710
711
712
713
714
715
716
717
718
719
720

31

741

a reduction in salinity, as quantified for March in Fig. 11e. Here, ρ anomalies averaged over the742

top 203m of the ocean in the Labrador Sea have been decomposed into the parts associated with743

salinity (ρS ) and temperature (ρT ) using an equation of state for sea water (McDougall et al. 2003),

744

and the salinity component is dominating. Associated with this is a reduction in convection in the

745

Labrador Sea, as depicted via the substantial reductions in March mixed layer depth in Fig. 11h,746

which uses the definition of Large et al. (1997).

747

We suspect, based on Fig. 11k that there is also an important role for positive salinity feedbacks

748

in the rapid AMOC decline with greenhouse gas forcing in CESM2. Figure 11k shows the lagged

749

regression of Labrador Sea density (ρ, ρS and ρT ) anomalies onto AMOC in the CESM2 pre-750

industrial control simulation using 10-year running means. We switch to using 10-year running751

means here for consistency with the study of Danabasoglu et al. (2019) which made use of this752

metric to indicate the role of density anomalies in driving and feeding back onto AMOC anomalies.

753

In the pre-industrial control variability, the maximum positive density anomalies due to changing

754

salinity lag the AMOC (dotted in Fig. 11k). This dominates over temperature feedbacks to lead755

to an overall positive density feedback on AMOC variability. We also suspect that differences in

756

this feedback between CESM2 and CESM1 are important in their differing AMOC behavior, to be

757

discussed in section 5.

758

Overall, when anthropogenic aerosol forcing is imposed on its own, it leads to an increase in sea

759

ice thickness (solid blue in Fig. 11b), a slight increase in Labrador Sea density through increased

760

salinity (solid blue in Fig. 11e) and a slight increase in AMOC strength (solid blue in Fig. 9e).

761

However, when the anthropogenic aerosol influence is inferred from LENS2−XAAER2, because762

the aerosols postpone the decline in sea ice thickness, ρ, mixed layer depth and AMOC in LENS2

763

compared to XAAER2, they lead to an apparent continued increase in Labrador Sea mixed layer

764

depth over the course of the simulation (Fig. 11h, blue dashed) and an increase in AMOC strength765(Fig. 9c, maroon dashed). The increase in AMOC strength is associated with enhanced ocean heat

766

transport into the NH high latitudes and, presumably in association with this, the sea ice thickness

767

starts to decline (blue dashed in Fig. 11b), as does the reduction in global mean temperature (Fig.

768

3b, maroon dashed).

32

(b) CESM2 sea ice thickness(a) AAER2 (LENS2 XAAER2)(c) CESM1 sea ice thickness

0.5

-2
-3

1.0

(kgm 3)(d) AAER2 (LENS2 XAAER2)

0.2

0.0

(kgm 3)0.2

Mixed layer depth (m)(g) AAER2 (LENS2 XAAER2)

400

200

0

200

Mixed Layer Depth (m)AMOC (Sv/Sv)

0.8

CESM2
CESM1

0.6
0.4
0.2
0.0
-0.2
-40 -30 -20 -10 0 10 20 30 40
Lag (years)

769
770
771
772
773
774
775
776
777
778
779
780
781

1980

2000

2020

-2
-3

2040

LENS1
XAAER1AAER1LENS1 XAAER1

1940

1960

1980

2000

2020

2040

0

0

-0.2

-0.2

-0.4

300
200
100
0
-100
-200
-300
-400
-500
-600

LENS2
LENS2 ( T)LENS2 ( S)

1940

1960

AAER2
AAER2 ( T)
AAER2 ( S)

1980

2000

2020

-0.8

2040

(h) CESM2 mixed layer depthLENS2
XAAER2AAER2LENS2 XAAER2

1940

1960

1980

2000

2020

0.05
S
0.04
T
0.03
0.02
0.01
0
-0.01
-0.02
leads
-0.03
-40 -30 -20 -10 0 10
Lag (years)

2040

30

300
200
100
0
-100
-200
-300
-400
-500
-600

LENS1
LENS1 ( T)LENS1 ( S)

1940

1960

AAER1
AAER1 ( T)
AAER1 ( S)1980

XAAER1XAAER1 ( T)
XAAER1 ( S)

2000

2020

2040

(i) CESM1 mixed layer depth

LENS1
XAAER1AAER1LENS1 XAAER1

1940

1960

1980

2000

2020

2040

(l) CESM1 Reg of onto AMOC

lags

20

-0.4
-0.6

XAAER2XAAER2 ( T)
XAAER2 ( S)(k) CESM2 Reg of onto AMOC(j) Reg of AMOC onto AMOC
(kgm 3Sv 1)

1.0

400

1960

-1

0.2 (f) CESM1 labrador sea density (top 203m)

-0.8

0.4

1940

0

0.2 (e) CESM2 labrador sea density (top 203m)

-0.6
0.4

LENS2
XAAER2AAER2LENS2 XAAER2(kgm 3)0.0-1Mixed layer depth (m)0.5

Sea ice thickness (m)0(kgm 3Sv 1)

1.0

1
Sea ice thickness (m)Sea ice thickness (m)

1

40

0.05
S
0.04
T
0.03
0.02
0.01
0
-0.01
-0.02
leads
-0.03
-40 -30 -20 -10 0 10
Lag (years)

lags

20

30

40

Fig. 11. (a) The difference in annual mean sea ice thickness anomalies between AAER2 and LENS2−XAAER2 during2030-2050. (b) 21-year running mean CESM2 annual mean sea ice thickness anomalies (relative to 1920-1940) averaged over220◦ E-360◦ E, 70◦ N-90◦ N (green region in (a)) for LENS2 (black), XAAER2 (pink), AAER2 (blue), and LENS2−XAAER2 (bluedashed), (c) as (b) but for CESM1. For AAER2 in (b) and LENS2−XAAER2 in (c) 95% confidence intervals are shown for a
3-member ensemble (light) and an ensemble of size equal to the one shown (dark) and for other experiments, the shading shows
the 95% confidence interval for a sample size equal to that in the ensemble. (d)-(f) show March density anomalies averaged over
the top 203 m of the ocean and the averaging region for (e) and (f) is the Labrador Sea (300◦ E-315◦ E, 53◦ N-65◦ N, green box in(d)). The LENS2−XAAER2 anomalies are not shown in (e) and (f) and instead the density anomalies in the other simulations aredecomposed into the part that is associated with temperature (ρT dashed) and the part that is associated with salinity (ρS dotted).
(g)-(i) are as (a)-(c) but for March mixed layer depth and the averaging in (h) and (i) is performed over the Labrador Sea. (j) lagged
auto-regression of annual mean AMOC. (k) lagged regression of March Labrador Sea density anomalies (and its temperature andsalinity components) in the top 203 m onto annual mean AMOC for CESM2. (l) as (k) but for CESM1. In (j)-(l) 10 year runningmeans are used for consistency with Danabasoglu et al. (2019)

33

782

c. Summary of method-dependence783

In this section we have investigated the dependence of the anthropogenic aerosol-forced [Ts ]784

response to the method used in CESM2 (AAER2 versus LENS2−XAAER2). The surface energy

785

balance indicates that aerosol forcing leads to a bigger decline in net TOA shortwave in AAER2

786

compared to LENS2−XAAER2 (Fig. 2c), which drives a bigger decline in global mean Ts and an

787

associated compensating decline in net TOA longwave radiation. The method-dependence of the788

TOA net shortwave radiation can be further narrowed down to difference in the surface upward789

shortwave radiation (Fig. 3g) linked to a difference in the surface shortwave albedo response (Fig.

790

3h).

791

We then discussed three potential sources of this albedo non-linearity. The first two are non-792

linearities or base state dependencies in the response of snow cover and sea ice fraction to Ts

793change. Snow cover over the continental regions surrounding the Arctic declines non-linearly794with warming. This is because the amount of time (in the summer at least) that is spent without795

any snow cover increases with warming and, as a result, the amount of time in which further

796

warming can influence the snow cover declines. The result is that the aerosol-forced cooling that797

is imposed within a cooler climate in AAER2, leads to a larger increase in snow cover than the798

aerosol-forced cooling that is imposed within a warmer climate that is influenced by CO2, as is the799

case in the LENS2 versus XAAER2 comparison. Sea ice exhibits a similar non-linear behavior

800

and there is also a dependence of the latitude at which sea ice grows with aerosol-forced cooling801

on the base state climate. In the cold climate within AAER2, additional sea ice grows at lower

802

latitudes and while the impacts of this on the global radiative balance is likely minimal during the803

summer, it may matter more during the spring and autumn when latitudinal gradients in insolation

804

are larger. These cryospheric effects lead to overall greater increases in albedo in AAER2 than805

in LENS2−XAAER2, and ultimately a differing global response as atmospheric and oceanic heat

806

transports respond.

807

Finally, there is also clearly a non-linear behavior of AMOC which further widens the discrep-

808

ancy between AAER2 and LENS2−XAAER2. In CESM2, with warming, the AMOC strength809

declines substantially and non-linearly. As a result, in LENS2 where aerosol forcing delays this810

AMOC decline, it leads to an apparently larger increase in AMOC strength due to aerosols in

811

LENS2−XAAER2 than is found due to aerosol forcing alone in AAER2. The associated increased

34

812

heat transport into the NH high latitudes in LENS2−XAAER2 would further alter sea ice and snow

813

cover with associated albedo changes.

814

5. Comparison between CESM2 and CESM1

815

We motivated the analysis in section 4 by aiming to determine the relative roles of the experimental

816

design used, the model physics and the aerosol forcing in leading to the differing anthropogenic

817

aerosol responses between CESM1 and CESM2. The results point to an important influence of the818

method in CESM2, but the same method-dependence was not found in CESM1. While the cooling819

due to anthropogenic aerosol forcing is greater in AAER1 than in LENS1−XAAER1 (solid versus

820

dashed teal in Fig. 3b), the difference is much smaller than in CESM2. There are two potential

821

reasons for this: (1) the model physics and dynamics in CESM1 is such that non-linearities are822

less important and (2) the fact that the imposed aerosol forcing declines more rapidly from the823

1980s (see the Appendix) may mean that CESM1 does not have as much of a chance for the824

non-linearities to lead to a big deviation between the methods. We cannot really explore the effect825

of (2) without simulations with CESM2 run under a lower aerosol emissions scenario and this826

would be a worthwhile avenue for future research to truly quantify the relative importance of these827

two factors. Nevertheless, differences between CESM1 and CESM2 do suggest that there is a role828

for model differences in leading to more non-linearity in CESM2 than in CESM1 and that it is not829

solely the difference in aerosol forcing that is responsible, as now discussed.

830

CESM1 exhibits much less of a difference in global upward surface shortwave between AAER1

831

and (LENS1−XAAER1) in the NH than was found between AAER2 and (LENS2−XAAER2),

832

although the SH difference is comparable (see supplemental Fig. 11, which is the equivalent of

833

Fig. 6 but for CESM1). CESM1 does not exhibit as substantial a method-dependence for sea ice

834

cover in the low latitudes of the Arctic or for snow cover over the continents surrounding the Arctic

835(supplemental Fig 12) but there are some similarities in the method-dependence for sea ice in the836

SH (supplemental Fig. 13). The differences between CESM1 and CESM2 in the NH can perhaps

837

be traced back to two differences between the CESM1 and CESM2 climates. In CESM1, the Arctic838

sea ice is thicker and more expansive (supplemental Fig. 14 and DuVivier et al. (2020)). A result of839

the more expansive sea ice in CESM1 is that there may be less room for sea ice to grow in the aerosol

840

only simulation before the continent is reached, limiting the differences in the extent to which the
35

841

sea ice fraction increases at lower latitudes in AAER1 compared to LENS1−XAAER1. The second

842

factor that may be important is that CESM1 has reduced summertime snow cover compared to843

CESM2 in portions of the continental regions surrounding the Arctic, in particular regions adjacent844

to Hudson Bay and to the south of the Barents-Kara Sea (supplemental Fig. 14c), which likely

845

makes those regions less non-linear in their snow cover response to temperature anomalies (recall

846

Fig. 7). These regions correspond to those where the wintertime snow density has increased the847

most in CESM2 compared to CESM1 in response to updated snow density parameterizations (see848

Fig. 5 of Simpson et al. (2022)), which is likely playing a role.

849

Another important difference between CESM1 and CESM2 is in the behavior of AMOC.

850

CESM1 does exhibit a method-dependence of the aerosol-forced AMOC response (Fig. 9d).

851

In LENS1−XAAER1 the AMOC strengthens more than in AAER1 (dashed versus solid in Fig.

852

9d) but this does not last for the full length of the simulation. The LENS1−XAAER1 AMOC

853

strength starts to decline substantially after about 1990 whereas the LENS2−XAAER2 AMOC

854

strengthening continues out to the end of the simulation (Fig. 9d versus c). This may be partly

855

due to the differing forcings between CESM1 and CESM2, but comparison of Figs. 9e and f856

makes clear that the behavior of AMOC in the XAAER simulations (pink), where aerosols are857

not evolving and greenhouse gas forcing is the primary driver, also differs considerably between

858

CESM1 and CESM2. In XAAER2, the AMOC declines much more rapidly after around 1980

859

than in XAAER1 and so the aerosol-forced strengthening inferred from the LENS-XAAER calcu-

860

lation is smaller during this period in CESM1 than in CESM2. Subsequently, when the aerosol

861

forcing starts to decline in the NH high latitudes, because it actually goes negative compared to862

the 1920-1940 baseline in CESM1, the aerosol forcing and GHG forcing act together to produce863

a sharper decline in AMOC in LENS1 than in XAAER1 in the 21st century (Fig. 9f black versus

864pink). The differing AMOC behavior is likely part of the reason why the high latitudes warm a lot

865

less in response to greenhouse gas forcing in GHG2 than in LENS1−XGHG1 (supplementary Fig.

86615). In CESM2, the greenhouse gas-forced decline in AMOC is greater than in CESM1, which867

reduces the northward heat transport into the high latitudes and reduces the warming there.

868

We speculate that an important factor in the differing AMOC responses between XAAER2 and869

XAAER1 is in the strength of salinity feedbacks. We can consider the AMOC decline to consist

870

of two parts: (1) the forcing which leads to the decline in the first place, which we argued for

36

871

CESM2 above was most likely the freshwater input to the Labrador Sea from sea ice melt; and, (2)872subsequent feedbacks which are triggered as the AMOC starts to decline, including the reduced873

advection of salty water from the low latitudes to the high latitudes. The various potential forcers

874

of AMOC decline (surface freshwater flux, surface heat flux and sea ice loss) can be compared875

between CESM1 and CESM2 for the XAAER experiment in supplemental Fig. 10. For both876

CESM1 and CESM2, sea ice loss appears as the most likely forcer of AMOC decline as it is877

the only one which leads the AMOC decline as opposed to lags it. However, a difference in878

sea ice loss cannot explain the differences in AMOC decline between XAAER1 and XAAER2,

879

because the sea ice loss is actually greater in XAAER1 than in XAAER2 (compare Figs. 11b880

and c), while the AMOC decline is greater in XAAER2. This suggests that the reason behind the881

difference in AMOC decline between XAAER1 and XAAER2 is more likely to be a difference882

in the feedbacks rather than in the initial forcing of AMOC decline. Figure 11j, which shows883

the lagged autoregression of AMOC onto itself within the pre-industrial control simulations of884

CESM1 and CESM2, demonstrates that the timescale of AMOC variability is longer in CESM2

885

than in CESM1. We may reasonably expect that a longer timescale AMOC variability is either due

886

to longer timescale forcing, whether that be through sea ice variability or surface flux variability,887

or due to stronger feedbacks onto AMOC variability which would lengthen the persistence of any

888

anomalies induced by the various forcers. Comparison of the lagged regression of density onto

889

AMOC in Figs. 11k and l shows that the salinity anomalies that lag AMOC in CESM1 are much

890

smaller than in CESM2, i.e., per unit Sverdrup increase in AMOC strength, the lagged increase891

in Labrador Sea salinity is greater in CESM2, which would provide a greater feedback onto an

892

AMOC change and, therefore, enhance the persistence of AMOC variability. This suggests that the893

positive salinity feedback onto AMOC anomalies may be stronger in CESM2 than in CESM1, for894

reasons that are currently unknown. This could lead to the more rapid AMOC decline in XAAER2,

895

even though the freshwater input through sea ice loss is smaller. A more detailed analysis of the896

AMOC decline in both simulations should be performed in future work to fully understand these897

differences. The updates to the ocean model in CESM2 compared to CESM1 are relatively minimal,898

but include the representation of mixing effects of estuaries, enhanced mesoscale eddy diffusivity

899

at depth, the use of prognostic chlorophyll for shortwave absorption, and the use of a salinity

900

dependent freezing point (Danabasoglu et al. 2020). Whether the differences in AMOC behavior

37

901

can be attributed to these ocean model changes or the changes to the coupled system introduced

902

through updates to the other components remains to be understood. Hassan et al. (2022) recently

903

argued that models that exhibit a greater AMOC response to forcing may do so because of a larger

904

feedback between the AMOC and cloud cover in the sub-polar North Atlantic. However, we find no905

evidence of a substantial feedback between the AMOC and sub-polar North Atlantic cloud cover in

906

CESM2 through regression of total cloud cover onto AMOC in the CESM2 pre-industrial control

907

simulation (not shown).

908

Overall, the comparison of the behavior of CESM1 and CESM2 makes clear that even though the909

single forcing experimental design matters within CESM2, it probably does so because of particular

910

features of both the representation of processes within the model and the imposed forcing.

911

6. Discussion and Conclusions

912

The implicit assumption when using single forcing experiments to attribute changes to individual

913

forcings is that non-linearities are negligible. As discussed in the introduction, prior studies have914

drawn mixed conclusions as to whether non-linearities are important. Some of the studies with

915

older model generations, e.g., Feichter et al. (2004) and Ming and Ramaswamy (2009), used a916

slab ocean, so any non-linearity related to AMOC would have been absent. The more recent

917

study of Deng et al. (2020) is the most relevant to the results presented here since they explored918

non-linearity within coupled CESM1 time-slice simulations. They did not find substantial non-919

linearities in GMST and TOA radiative fluxes, aligned with our findings that the “only” versus

920“all-but-one” method does not dramatically alter those responses to aerosol forcing in transient921

experiments with CESM1. Further probing of other features by Deng et al. (2020) did reveal other922

non-linearities, specifically in September-November Arctic sea ice decline and in summertime923

precipitation over East Asia. The sense of their sea ice non-linearity was that when greenhouse924

gases and aerosols were imposed together, there was less sea ice decline than summing up the925

contributions from greenhouse gases and aerosols separately, which is the opposite of what we find926

for CESM2.

927

The considerable non-linearities that we infer in CESM2 from the difference between the “only”

928

and “all-but-one” methods for anthropogenic aerosols and comparison to the behavior in CESM1

929

and these other previous studies makes clear that non-linearities in the response to forcings can
38

930

be highly dependent on the model physics and/or the forcings used. Indeed Menary et al. (2020)

931

find that in CMIP6 models in general, the aerosol and greenhouse gas-forced AMOC anomalies932

do approximately sum up to the response when all forcings are applied together. CESM2 is clearly933

a more non-linear model than CESM1, particularly when it comes to the AMOC response to934

forcings, but also likely in the impact of surface shortwave albedo feedbacks. For the sea ice

935

aspects, the version of CESM2 used here is known to be deficient in its representation of sea ice

936(DuVivier et al. 2020) so sea ice changes should be interpretted with caution. For snow cover,937

further investigation is required to determine whether summertime snow cover in CESM2 is more

938aligned with observations than CESM1 (supplemental Fig. 14a-c), although Wieder et al. (2022)

939

indicate that CESM2 does have too much snow water equivalent in the springtime in the regions

940

adjacent to Hudson Bay and to the south of the Kara Sea. Much work also remains to be done to941

fully understand the differences in AMOC variability and change between CESM1 and CESM2

942

and to determine whether we trust one more than the other.

943

Overall, the method-dependence found for the aerosol-forced response in CESM2, raises the944question, what is the more appropriate method to use in single forcing experiments? Our experience945

with two generations of CESM indicates that the method used may matter for some models and/or

946

forcings more than others. Ultimately, there is probably no getting around the non-linearities that947

exist in CESM2, particularly that due to the AMOC and one method is not going to necessarily

948

give you a more correct answer than the other. Using the “only” method we would conclude that949

greenhouse gases are giving rise to a dramatic decline in AMOC strength that starts in the mid-20th

950

century, whereas the reality is that the greenhouse gases don’t have this same effect when they951

are imposed together with aerosol forcing. Using the “all-but-one” method we would conclude952

that aerosols give rise to an increasing AMOC strength at least out to 2050, but the reality is that953

they are only apparently doing that because they have prevented the greenhouse gas-forced AMOC

954

decline. The AMOC response to forcings is non-linear and neither method alone would provide955

the complete picture and we should be aware of such non-linearities in our interpretation.

956

This new CESM2 dataset has been released to the research community and we expect that there957

are many more interesting insights that can be gained from it. We also expect that further insights

958

can be gained by building on the dataset provided here through modified experimental design959

and/or forcing combinations. Unlike the CESM1 single forcing large ensemble, this new ensemble

39

960

offers the opportunity to assess additivity of the different forcing contributions in comparison to the961

overall LENS2 response. The results presented here highlight the importance of non-linearities in

962

interpreting single-forcing simulations, while simultaneously highlighting pertinent mechanisms

963

underlying these non-linearities that may be of value for future endeavors. It is our hope that future964

work will make use of this dataset to further explain the role of individual forcings and identify

965

their interactions in the evolution of the Earth system.

40

966

Acknowledgments. This work is supported by the National Center for Atmospheric Research,967

which is a major facility sponsored by the National Science Foundation under the Cooperative

968

Agreement 1852977. Computing and data storage resources, including the Cheyenne supercom-

969

puter (https://doi.org/10.5065/D6RX99HX), were provided by the Computational and Information

970

Systems Laboratory at NCAR. We also acknowledge the Community Earth System Model version

971

1 Large Ensemble Community Project, the Community Earth System Model version 2 Large En-

972

semble Community Project and supercomputing resources provided by the IBS Center for Climate

973

Physics in South Korea (https://doi.org/10.5194/esd-2021-50). The CESM2 large ensemble simu-

974

lations were conducted on the IBS/ICCP supercomputer “Aleph”, a 1.43 petaflop high-performance

975

Cray XC50-LC Skylake computing system. The work of KBR was supported by the Institute for976

Basic Sciences (IBS), Republic of Korea, under IBS-R028-D1. IRS is grateful to Brian Medeiros,

977

John Fasullo, Robb Jnglin Wills, Andrew Gettelman and Marika Holland for helpful discussions.

978

Data availability statement.

979

able through NCAR’s Climate Data Gateway with guidance for data access located980

at the following website https://www.cesm.ucar.edu/projects/community-projects/LENS/data-981

sets.html.

982

mate Data Gateway with guidance for data access located at the following website

983

https://www.cesm.ucar.edu/projects/community-projects/LENS2/data-sets.html.

984

single forcing large is available through the Climate Data Gateway as described at985

https://www.cesm.ucar.edu/working-groups/climate/simulations/cesm2-single-forcing-le.

986

analysis987

https://github.com/islasimpson/singleforcing_paper and the data required to make the figures are

988

available at https://doi.org/10.5065/ayy9-fr14

The CESM1 datasets used in this study are freely avail-

The CESM2 large ensemble is also available through NCAR’s Cli-

codes

for

reproducing

thefiguresinthismanuscriptare

The CESM2

available

The
at

989

APPENDIX

990

Comparison of the anthropogenic aerosol forcing between the CESM2 and CESM1 single

991

forcing large ensembles

992

Given that much of the analysis in this study focusses on the anthropogenic aerosol-forced993

response, we provide a comparison of the aerosol forcing between AAER2 and AAER1 in Fig. A1

994

with a focus on two species (BC and SO4 ), while the emissions and burdens of other species can be
41

995

found in supplemental Figs. 3 and 4. A comparison of Figs. A1a and A1b reveals that the trends in996

anthropogenic aerosol optical depth at 550 nm (AOD) over 2000-2050 are very different between

997

the AAER2 (historical to SSP3-7.0) and AAER1 (historical to RCP8.5) simulations. AAER1 shows

998

declines in AOD over eastern North America, Europe and China with relatively small increases

999

in AOD over Africa and India (Fig. A1b). In AAER2, the declines over eastern North America

1000

and Europe are much smaller compared to AAER1 and the AOD increases over China in AAER2

1001

while it decreases in AAER1. Over Africa and India, AOD is increasing much more in AAER2

1002

than in AAER1. This difference primarily stems from the difference in emissions, but differences1003

in the model physics also play a role in the differing overall aerosol burdens between AAER1 and1004

AAER2, as now discussed.

1005

The global BC emissions in AAER1 and AAER2 are fairly similar until about the year 2000,

1006

but after that they increase in AAER2 and decrease in AAER1 (Fig. A1 c) (a similar trajectory is1007

seen for POM in supplemental Fig. 3d). SO4 emissions are also similar until about the year 2000,

1008

but then they remain fairly steady in AAER2 while declining in AAER1 (Fig. A1d) (a similar

1009

trajectory is seen for SO2 in supplemental Fig. 3c). For BC, while the emissions are slightly lower

1010

over the 19th and early 20th centuries in AAER2 compared to AAER1 (Fig. A1c), the BC burden

1011

is higher in AAER2 compared to AAER1 (Fig. A1e). The lifetime of BC, estimated by the ratio

1012

of the global BC burden to the global BC deposition flux, is longer in CESM2 (6.41 days) than in

1013

CESM1 (3.64 days, Fig. A1g). This nearly two-fold increase in BC lifetime can be understood as

1014

resulting from the BC wet deposition flux associated with a given burden being smaller in CESM2

1015

than in CESM1 and the wet deposition flux has changed because the representation of the aging

1016

of primary carbonaceous aerosols in CESM2 delays BC removal via wet deposition (section 2.a.1

1017

and see the differing wet deposition rates at a given global burden in Fig. A1g). For SO4 , the1018

emissions are rather similar over the historical period (Fig. A1d) but the burden (Fig. A1f) is

1019

higher in AAER1 than in AAER2 for reasons that are not totally clear given that the deposition1020

rates are fairly comparable between CESM1 and CESM2 (Fig. A1h).

1021

In summary, the difference in BC and SO4 emissions is the primary contributor to the difference1022

in burden (and associated AOD) trends between AAER2 and AAER1, with some additional

1023

modification due to the differing model physics. The SSP3-7.0 scenario is a higher aerosol

1024

emission scenario (Gidden et al. 2019) than the CMIP5 RCP8.5 scenario used in the CESM1 single

42

1025

forcing large ensemble and the emissions also differ slightly over the historical period as emissions

1026

inventories were revised between CMIP5 and CMIP6, although it should be noted that it has been1027

argued that the increasing emissions over eastern China in the last decade of the historical period in1028

the CMIP6 emissions are incorrect (Wang et al. 2021). The impact of such differences in emissions1029

can end up being as large as the impact of changing from one model version to the next (e.g., Fyfe

1030

et al. 2021).

43

1037

References

1038

Allen, R. J., J. R. Norris, and M. Kovilakam, 2014: Influence of anthropogenic aerosols and the1039

Pacific Decadal Oscillation on tropical belt width. Nature Geoscience, 7, 270–274.

1040

Baek, S. H., Y. Kushnir, M. Ting, J. E. Smerdon, and J. M. Lora, 2022: Regional Signatures of1041

Forced North Atlantic SST Variability: A Limited Role for Aerosols and Greenhouse Gases.

1042

Geophysical Research Letters, 49, https://doi.org/https://doi.org/10.1029/2022GL097794.

1043

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

2018: The path to CAM6: coupled simulations with CAM5.4 and CAM5.5. Geosci. Model

1045

Dev., 11, 235–255.

1046

Bonfils, C. J. W., B. D. Santer, J. C. Fyfe, K. Marvel, T. J. Phillips, and S. R. H. Zimmerman,

1047

2020: Human influence on joint changes in temperature, rainfall and continental aridity. Nature

1048

Climate Change, 10, 726–731.

1049

Brodeau, L., and T. Koenigk, 2016: Extinction of the northern oceanic deep convection in an1050

ensemble of climate model simulations of the 20th and 21st centuries. Clim. Dyn., 46, 2863–

1051

2882.

1052

Capotondi, A., C. Deser, A. S. Phillips, Y. Okumura, and S. M. Larson, 2020: ENSO and Pacific1053

decadal variability in the Community Earth System Model version 2. J. Adv. Mod. Earth Sys.,

1054

12, e2019MS002 022.

1055

Chiang, F., O. Mazdiyasni, and A. AghaKouchak, 2021: Evidence of anthropogenic impacts

1056

on global drought frequency, duration, and intensity. Nature Communications, 12, 2754,

1057

https://doi.org/https://doi.org/10.1038/s41467-021-22314-w.

1058

Dagan, G., P. Stier, and D. Watson-Parris, 2020: Aerosol Forcing Masks and Delays the Formation1059

of the North Atlantic Warming Hole by Three Decades. Geophys. Res. Lett., 47, https://doi.org/

1060

https://doi.org/10.1029/2020GL090778.

1061

1062

Danabasoglu, G., S. C. Bates, B. P. Briegleb, S. R. Jayne, and W. G. Jochum, M. amd Large, 2012:
The CCSM4 ocean component. J. Clim., 25, 1361–1389.
44

(a) AOD trend 2000-2050, AAER2

0.08

0.06

0.04

(b) AOD trend 2000-2050, AAER1

0.02

0.00

4
2

1033
1034
1035
1036

Burden (×10 6 kgm 2)(g) Black Carbon

3.64 days
LENS1
LENS2
LENS1 (dry)LENS1 (wet)LENS2 (dry)LENS2 (wet)06.41 days

0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0
Burden (×10 7 kgm 2)(f) SO4 Burden

3.0

AAER1
AAER2

2.5
2.0
1.5
1.0
0.5

Deposition (×10 12kgm 2s 1)

Burden (×10 7 kgm 2)(e) Black Carbon Burden
4.0
AAER1
3.5
AAER2
3.0
2.5
2.0
1.5
1.0
0.5
0
1850 1875 1900 1925 1950 1975 2000 2025 2050

Deposition (×10 13kgm 2s 1)
1032

0.08

AAER1
6
AAER2
5
4
3
2
1
0
1850 1875 1900 1925 1950 1975 2000 2025 2050

1850 1875 1900 1925 1950 1975 2000 2025 2050

1031

0.06

Emissions (Tg/y)Emissions (Tg/y)

6

7
6
5
4
3
2
1
0

0.04

(d) SO4 emissions(c) Black Carbon emissions
AAER1
AAER2

8

0.02

AOD trends (per 51y)

0
12
11
10
9
8
7
6
5
4
3
2
1
0

1850 1875 1900 1925 1950 1975 2000 2025 2050(h) SO4
LENS1
LENS2
LENS1 (dry)LENS1 (wet)LENS2 (dry)LENS2 (wet)00.5

3.37 days

1.0 1.5 2.0 2.5
Burden (×10 6 kgm 2)

3.71 days

3.0

3.5

Fig. A1. (a) and (b) Trends in annual mean AOD at 550nm between 2000 and 2050 of AAER2 and AAER1, respectively. (c)Global annual mean emissions of BC for AAER1 (teal) and AAER2 (maroon). (d) same as (c) but for SO4 aerosol. (e) Global
annual mean BC burden for AAER1 (teal) and AAER2 (maroon). (f) as (e) but for SO4 . (g) Global annual mean deposition fluxes
versus burden for a single member of LENS2 and LENS1 (used rather than the single forcing experiments because the deposition
fluxes were not output in the CESM1 single forcing large ensemble). Circles show the full deposition flux (dry + wet), stars show
the wet deposition flux and diamonds show the dry deposition flux. (h) as (g) but for SO4 .

45

1063

Danabasoglu, G., L. Landrum, S. G. Yeager, and P. R. Gent, 2019: Robust and Nonrobust Aspects

1064

of Atlantic Meridional Overturning Circulation Variability and Mechanisms in the Community

1065

Earth System Model. J. Clim., 32, 7349–7368.

1066

1067

1068

1069

Danabasoglu, G., and Coauthors, 2020: The Community Earth System Model 2 (CESM2). J. Adv.
Mod. Earth Sys., 12, https://doi.org/10.1029/2019MS001916.
Delworth, T. L., and F. Zeng, 2014: Regional rainfall decline in Australia attributed to anthropogenic greenhouse gases and ozone levels. Nature Geoscience, 7, 583–587.

1070

Delworth, T. L., F. Zeng, L. Zhang, R. Zhang, G. Vecchi, and X. Yang, 2017: The Central Role of

1071

Ocean Dynamics in Connecting the North Atlantic Oscillation to the Extratropical Component

1072

of the Atlantic Multidecadal Oscillation. J. Clim., 30, 3789–3805.

1073

1074

Deng, J., A. Dai, and H. Xu, 2020: Nonlinear Climate Response to Increasing CO2 and Anthropogenic Aerosols Simulated by CESM1. J. Clim., 33, 281–301.

1075

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

1076

due to CMIP6 biomass burning emissions. Science Advances, 8 (30), https://doi.org/10.1126/

1077

sciadv.abo2405.

1078

1079

1080

1081

Deser, C., R. A. Tomas, and L. Sun, 2015: The Role of Ocean-Atmosphere Coupling in the Zonal
Mean Atmospheric Response to Arctic Sea Ice Loss. J. Clim., 28, 2168–2186.
Deser, C., and Coauthors, 2020a: Insights from Earth system model initial-condition large ensembles and future prospects. Nat. Clim. Change, 10, 277–286.

1082

Deser, C., and Coauthors, 2020b: Isolating the Evolving Contributions of Anthropogenic Aerosols

1083

and Greenhouse Gases: A New CESM1 Large Ensemble Community Resource. J. Clim., 33,1084

7835–7858.

1085

Dittus, A. J., E. Hawkins, J. I. Robson, D. M. Smith, and J. Wilcox, L, 2021: Drivers of Recent

1086

North Pacific Decadal Variability: The Role of Aerosol Forcing. Earth’s Future, 9, https://doi.org/

1087

https://doi.org/10.1029/2021EF002249.
46

1088

Dixon, K. W., T. L. Delworth, M. J. Spelman, and R. J. Stouffer, 1999: The influence of transient1089

surface fluxes on North Atlantic overturning in a coupled GCM Climate Change Experiment.

1090

Geophys. Res. Lett., 26, 2749–2752.

1091

1092

Dong, B., R. T. Sutton, E. Highwood, and L. Wilcox, 2014: The Impacts of European and Asian
Anthropogenic Sulfur Dioxide Emissions on Sahel Rainfall. J. Clim., 27, 7000–7017.

1093

DuVivier, A. K., M. M. Holland, J. E. Kay, S. Tilmes, A. Gettelman, and D. A. Bailey, 2020:

1094

Arctic and Antarctic sea ice mean state in the Community Earth System Model version 2 and1095

the influence of atmospheric chemistry. J. Geophys. Res.: Oceans, 125, e2019JC015 934.

1096

England, M. R., I. Eisenman, N. J. Lutsko, and T. J. W. Wagner, 2021: The Recent Emer-

1097

gence of Arctic Amplification. Geophys. Res. Lett., 48, https://doi.org/https://doi.org/10.1029/

1098

2021GL094086.

1099

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

1100

Overview of the Coupled Model Intercomparison Project Phase 6 (CMIP6) experimental design1101

and organization. Geoscientific Model Development, 9, 1937–1958.

1102

Fasullo, J. T., P. R. Gent, and R. S. Nerem, 2020: Sea Level Rise in the CESM Large Ensemble:

1103

The Role of Individual Climate Forcings and Consequences for the Coming Decades. J. Clim.,

1104

33, 6911–6927.

1105

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

1106

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

1107

Biomass Burning Emissions. Geophys. Res. Lett., 49, e2021GL097 420, https://doi.org/https:

1108//doi.org/10.1029/2021GL097420.

1109

1110

Feichter, J., E. Roeckner, U. Lohmann, and B. Liepert, 2004: Nonlinear Aspects of the Climate
Response to Greenhouse Gas and Aerosol Forcing. J. Clim., 17, 2384–2398.

1111

Fyfe, J. C., V. V. Kharin, B. D. Santer, J. N. S. Cole, and N. P. Gillett, 2021: Significant impact1112

of forcing uncertainty in a large ensemble of climate model simulations. Proc. Natl. Acad. Sci.

1113

(USA), 118, e2016549 118.

1114

1115

Gettelman, A., and H. Morrison, 2015: Advanced two-moment bulk microphysics for global
models. Part 2: Off-line tests and comparison with other schemes. J. Clim, 28, 1268–1287.
47

1116

1117

1118

1119

Gettelman, A., and Coauthors, 2019: The Whole Atmosphere Community Climate Model Version
6 (WACCM6). J. Geophys. Res.: Atm, 124, https://doi.org/10.1029/2019JD030943.
Giannini, A., and A. Kaplan, 2019: The role of aerosols and greenhouse gases in Sahel drought
and recovery. Climatic Change, 152, 449–466.

1120

Gidden, M. J., and Coauthors, 2019: Global emissions pathways under different socioeconomic

1121

scenarios for use in CMIP6: a dataset of harmonized emissions trajectories through the end of1122

the century. Geosci. Model Dev., 12, 1443–1475.

1123

1124

1125

1126

Gillett, N. P., and Coauthors, 2016: The Detection and Attribution Model Intercomparison Project(DAMIP v1.0) contribution to CMIP6. Geosci. Model Dev., 9, 3685–3697.
Gillett, N. P., and Coauthors, 2021: Constraining human contributions to observed warming since
the pre-industrial period. Nat. Clim. Change, 11, 207–212.

1127

Hassan, T., R. J. Allen, W. Liu, and C. A. Randles, 2021: Anthropogenic aerosol forcing of the1128

Atlantic meridional overturning circulation and the associated mechanisms in CMIP6 models.

1129

Atmos. Chem. Phys., 21, 5821–5846.

1130

Hassan, T., and Coauthors, 2022: Air quality improvements are projected to weaken the Atlantic

1131

meridional overturning circulation through radiative forcing effects. Nature Communications

1132

Earth and Environment, 3, https://doi.org/10.1038/s43247-022-00476-9.

1133

Hirasawa, H., P. J. Kushner, M. Sigmond, J. Fyfe, and C. Deser, 2020: Anthropogenic Aerosols

1134

Dominate Forced Multidecadal Sahel Precipitation Change through Distinct Atmospheric and1135

Oceanic Drivers. J. Clim., 33, 10 187–10 204.

1136

Hunke, E. C., W. H. Lipscomb, A. K. Turner, N. Jeffery, and S. Elliott, 2015: Cice: The los alamos

1137

sea ice model. documentation and software user’s manual. version 5.1. Tech. Rep. LA-CC-06-

1138

012, iT-3 Fluid Dynamics Group, Los Alamos National Laboratory.

1139

1140

Hurrell, J. W., and Coauthors, 2013: The Community Earth System Model: A Framework for
Collaborative Research. Bull. Amer. Meteor. Soc., 94, 1339–1360.

1141

IPCC, 2021: Summary for Policy Mackers. In: Climate Change 2021: The Physical Science Basis.

1142

Contribution of Working Group 1 to the Sixth Assessment Report of the Intergovernmental Panel
48

1143

on Climate Change. Tech. rep., Cambridge University Press, Cambridge, United Kingdom and1144

New York, NY, USA, 3-32 pp. https://doi.org/10.1017/9781009157896.001.

1145

Jahn, A., and M. M. Holland, 2013: Implications of Arctic sea ice changes for North Atlantic deep

1146

convection and the meridional overturning circulation in CCSM4-CMIP5 simulations. Geophys.

1147

Res. Lett., 40, 1206–1211.

1148

1149

Kay, J. E., and Coauthors, 2015: The Community Earth System Model (CESM) Large Ensemble
Project. Bull. Amer. Meteor. Soc., 96, 1333–1349.

1150

Lamarque, J. F., G. P. Kyle, M. Meinshausen, K. Riahi, S. J. Smith, D. P. van Vuuren, A. J. Conley,1151

and F. Vitt, 2011: Global and regional evolution of short-lived radiatively active gases and1152

aerosols in the Representative Concentration Pathways. Climatic Change, 109, 191–212.

1153

Large, W. G., G. Danabasoglu, S. C. Doney, and J. C. McWilliams, 1997: Sensitivity to surface

1154

forcing and boundary layer mixing in the NCAR CSM ocean model: Annual-mean climatology.

1155

J. Phys. Oceanogr., 27, 2418–2447.

1156

Lawrence, D., and Coauthors, 2019: The Community Land Model, version 5: Description of1157

new features, benchmarking and impact of forcing uncertainty. J. Adv. Mod. Earth Sys., 11,1158

4245–4287.

1159

Lawrence, D. M., and Coauthors, 2011: Parameterization improvements and functional and1160

structural advances in Version 4 of the Community Land Model. J. Adv. Mod. Earth Sys.,

1161

3, https://doi.org/https://doi.org/10.1029/2011MS00045.

1162

1163

Li, H., A. Federov, and W. Liu, 2021: AMOC Stability and Diverging Response to Arctic Sea Ice
Decline in Two Climate Models. J. Clim., 34, 5443–5460.

1164

Li, X., M. Ting, and D. E. Lee, 2018: Fast Adjustments of the Asian Summer Monsoon to1165

Anthropogenic Aerosols. Geophys. Res. Lett., 45, 1001–1010, https://doi.org/https://doi.org/10.

1166

1002/2017GL076667.

1167

Liu, X., P.-L. Ma, H. Wang, S. Tilmes, B. Singh, R. C. Easter, S. J. Ghan, and P. J. Rasch, 2016:

1168

Description and evaluation of ta new four-mode version of the Model Aerosol Module (MAM4)1169

within version 5.3 of the Community Atmosphere Model. Geosci. Model Dev., 9, 505–522.
49

1170

Liu, X., and Coauthors, 2012: Toward a minimal representation of aerosols in climate models:1171

description and evaluation in the Community Atmosphere Model CAM5. Geosci. Model Dev.,

1172

5, 709–739.

1173

Lombardozzi, D. L., Y. Lu, P. J. Lawrence, D. M. Lawrence, S. Swenson, K. W. Oleson, W. R.

1174

Wieder, and E. A. Ainsworth, 2020: Simulating agriculture in the community land model version

1175

5. J. Geophys. Res.: Biogeosciences, 125, e2019JG005 529.

1176

1177

Manabe, S., and R. J. Stouffer, 1993: Century-scale effects of increased atmospheric CO2 on theocean-atmosphere system. Nature, 364, 215–217.

1178

McDougall, T. J., D. R. Jackett, D. G. Wright, and R. Feistel, 2003: Accurate and Computationally1179

Efficient Algorithms for Potential Temperature and Density of Seawater. Journal of Atmospheric

1180

and Oceanic Technology, 20, 730–741.

1181

Meehl, G. A., C. Shields, J. M. Arblaster, H. Annamalai, and R. Neale, 2020: Intraseasonal,

1182

Seasonal, and Interannual Characteristics of Regional Monsoon Simulations in CESM2. J. Adv.

1183

Mod. Earth Sys., 12, https://doi.org/https://doi.org/10.1029/2019MS001962.

1184

Meehl, G. A., W. M. Washington, C. M. Ammann, J. M. Arblaster, T. M. L. Wigley, and C. Tebaldi,1185

2004: Combinations of Natural and Anthropogenic Forcings in Twentieth-Century Climate. J.

1186

Clim., 17, 3721–3727.

1187

1188

1189

1190

Meinshausen, M., and Coauthors, 2011: The RCP greenhouse gas concentrations and their extensions from 1765 to 2300. Climatic Change, 109, 213–241.
Meinshausen, M., and Coauthors, 2020: The shared socio-economic pathway (SSP) greenhousegas concentrations and their extensions to 2500. Geosci. Model Dev., 13, 3571–3605.

1191

Menary, M. B., J. Robson, R. P. Allan, B. B. B. Booth, C. Cassou, and G. Gastineau, 2020:

1192

Aerosol-forced AMOC changes in CMIP6 historical simulations. Geophys. Res. Lett., 47,1193

e2020GL088 166.

1194

1195

Ming, Y., and V. Ramaswamy, 2009: Nonlinear Climate and Hydrological Responses to Aerosol
Effects. J. Clim., 22, 1329–1339.
50

1196

Monerie, P.-A., L. J. Wilcox, and A. G. Turner, 2022: Effects of Anthropogenic Aerosol and

1197

Greenhouse Gas Emissions on Northern Hemisphere Monsoon Precipitation: Mechanisms and

1198

Uncertainty. J. Clim., 35, 2305–2326.

1199

Morrison, H., and A. Gettelman, 2008: A new two-moment bulk stratiform cloud microphysics

1200

scheme in the NCAR Community Atmosphere Model (CAM3). Part 1: Description and numer-

1201

ical tests. J. Clim., 21, 3642–3659.

1202

Mueller, B. L., N. P. Gillett, A. H. Monahan, and F. W. Zwiers, 2018: Attribution of Arctic Sea

1203

Ice Decline from 1953 to 2012 to Influences from Natural, Greenhouse Gas and Anthropogenic

1204

Aerosol Forcing. J. Clim., 31, 7771–7787.

1205

Pendergrass, A. G., D. B. Coleman, C. Deser, F. Lehner, N. Rosenbloom, and I. R. Simpson, 2019:

1206

Nonlinear Response of Extreme Precipitation to Warming in CESM1. Geophys. Res. Lett., 46,1207

10 551–10 560.

1208

Robson, J., and Coauthors, 2022: The role of anthropogenic aerosol forcing in the 1850-1985

1209

strengthening of the AMOC in CMIP6 historical simulations. J. Clim., https://doi.org/https:

1210//doi.org/10.1175/JCLI-D-22-0124.1.

1211

1212

Rodgers, K. B., and Coauthors, 2021: Ubiquity of human-induced changes in climate variability.
Earth Sys. Dynam., 12, 1393–1411.

1213

Seong, M., S. Min, Y. Kim, X. Zhang, and Y. Sun, 2021: Anthropogenic Greenhouse Gas and1214

Aerosol Contributions to Extreme Temperature Changes during 1951-2015. J. Clim., 34, 857–

1215

870.

1216

1217

Shi, J.-R., Y.-O. Kwon, and S. Wijffels, 2022: Two Distinct Modes of Climate Responses to the
Anthropogenic Aerosol Forcing Changes. J. Clim., 35, 3445–3457.

1218

Simpson, I. R., D. M. Lawrence, S. C. Swenson, C. Hannay, K. A. McKinnon, and J. E. Truesdale,1219

2022: Improvements in Wintertime Surface Temperature Variability in the Community Earth

1220

System Model Version 2 (CESM2) Related to the Representation of Snow Density. J. Adv. Mod.

1221

Earth Sys., 14, e2021MS002 880.
51

1222

Simpson, I. R., and Coauthors, 2020: An evaluation of the large-scale atmospheric circulation

1223

and its variability in CESM2 and other CMIP models. J. Geophys. Res.: Atmospheres, 125,1224

e2020JD032 835.

1225

Singh, D., S. P. McDermid, B. I. Cok, M. J. Puma, L. Nazarenko, and M. Kelley, 2018: Distinct

1226

Influences of Land Cover and Land Management on Seasonal Climate. J. Geophys. Res. Atm.,

1227

123, 12 017–12 039, https://doi.org/https://doi.org/10.1029/2018JD028874.

1228

Smith, D., and Coauthors, 2022: Attribution of multi-annual to decadal changes in the climate sys-1229

tem: The Large Ensemble Single Forcing Model Intercomparison Project (LESFMIP). Frontiers

1230

in Climate, https://doi.org/10.3389/fclim.2022.955414.

1231

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

1232

component of the Community Climate System Model (CCSM). Tech. rep., LANL Tech. Report.,

1233

141 pp. https://doi.org/LAURâĂŘ10âĂŘ01853.

1234

1235

Tilmes, S., and Coauthors, 2019: Climate Forcing and Trends of Organic Aerosols in the Community Earth System Model (CESM2). J. Adv. Mod. Earth Sys., 11, 4323–4351.

1236

Touma, D., S. Stevenson, F. Lehner, and S. Coats, 2021: Human-driven greenhouse gas and aerosol

1237

emissions cause distinct regional impacts on extreme fire weather. Nature Communications, 12,1238

https://doi.org/https://doi.org/10.1038/s41467-020-20570-w.

1239

Undorf, S., D. Polson, M. A. Bollasina, Y. Ming, A. Schurer, and G. C. Hegerl, 2018: Detectable

1240

Impact of Local and Remote Antrhopogenic Aerosols on the 20th Century Chnges of West1241

African and South Asian Monsoon Precipitation. J. Geophys. Res. Atm., 123, 4871–4889.

1242

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

1243

merging satellite observations with proxies and fire models (1750-2015). Geosci. Model Dev.

1244

Discuss., 10, 3329–3357, https://doi.org/https://doi.org/10.5194/gmd-10-3329-2017.

1245

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

1246

Loss of Arctic Sea Ice via Ocean Dynamics. Geophys. Res. Lett., 45, https://doi.org/https:

1247//doi.org/10.1029/2018GL077325.

1248

1249

Wang, Z., and Coauthors, 2021: Incorrect Asian aerosols affecting the attribution and projectionof regional climate change in CMIP6 models. npj Clim Atmos Sci, 4.
52

1250

1251

Watanabe, M., and H. Tatebe, 2019: Reconciling roles of sulphate aerosol forcing and internalvariability in Atlantic multidecadal climate changes. Climate Dynamics, 53, 4651–4665.

1252

Weaver, A. J., M. Eby, M. Kienast, and O. A. Saenko, 2007: Response of the Atlantic merid-1253

ional overturning circulation to increasing atmospheric CO2: sensitivity fo mean climate state.

1254

Geophys. Res. Lett., 34, L05708.

1255

Wieder, W. R., D. Kennedy, F. Lehner, K. N. Musselman, K. B. Rodgers, N. Rosenbloom, I. R.

1256

Simpson, and R. Yamaguchi, 2022: Pervasive alterations to snow-dominated ecosystem func-

1257

tions under climate change. Proc. Nat. Ac. Sci, 119, https://doi.org/https://doi.org/10.1073/pnas.

1258

2202393119.

1259

Yang, Q., T. H. Dixon, P. G. Myers, J. Bonin, D. Chambers, M. R. van den Broeke, M. H.

1260

Ribergaart, and J. Mortensen, 2016: Recent increases in Arctic freshwater flux affects Labrador

1261

Sea convection and Atlantic overturning circulation. Nature Communications, 7, 10 525.

1262

Zhang, S., P. Stier, G. Dagan, and M. Wang, 2021: Anthropogenic Aerosols Modulated 20th-

1263

Century Sahel Rainfall Variability Via Their Impacts on North Atlantic Sea Surface Temperature.

1264

Geophys. Res. Lett., 49, https://doi.org/https://doi.org/10.1029/2021GL095629.

53