Using the fast impact of anthropogenic aerosols on regional land temperature to constrain aerosol forcing

See allHide authors and affiliations

Science Advances  05 Aug 2020:
Vol. 6, no. 32, eabb5297
DOI: 10.1126/sciadv.abb5297


Anthropogenic aerosols have been postulated to have a cooling effect on climate, but its magnitude remains uncertain. Using atmospheric general circulation model simulations, we separate the land temperature response into a fast response to radiative forcings and a slow response to changing oceanic conditions and find that the former accounts for about one fifth of the observed warming of the Northern Hemisphere land during summer and autumn since the 1960s. While small, this fast response can be constrained by observations. Spatially varying aerosol effects can be detected on the regional scale, specifically warming over Europe and cooling over Asia. These results provide empirical evidence for the important role of aerosols in setting regional land temperature trends and point to an emergent constraint that suggests strong global aerosol forcing and high transient climate response.


The global land surface air temperature (LSAT) has increased markedly since the 1960s (1, 2), and there is growing evidence that anthropogenic forcings have made a substantial contribution to this warming (3). The total response of land temperature to forcing (F) can be decomposed into two componentsδLSAT=(δLSATδFO+δLSATδOFδOδF)δF(1)

The first term on the right-hand side of Eq. 1 represents the atmospheric/land-only, fast response, and the second term represents the ocean (O)–mediated, slow response. Both responses are present in coupled atmosphere-ocean general circulation models (AOGCMs), the primary tool for detecting and attributing the influence of anthropogenic forcings on surface temperature (46). The slow response contains much of the uncertainty in AOGCM-simulated climate sensitivity, and it is this component that is most affected by low-frequency natural variability, thus limiting the robustness of AOGCM results (6, 7). In contrast, the fast response, which contains valuable information on forcing, is less affected by these complications.

Atmospheric general circulation models (AGCMs) forced by observed oceanic boundary conditions (sea surface temperatures and sea ice concentrations) (8) can be used to separate the fast from the slow response by circumventing the issues pertaining to climate sensitivity and natural variability. AGCMs have been used to understand the changes in the atmospheric circulation and hydrological cycle (9, 10) but have received little attention in the analysis of land temperature, perhaps due to the concern that the oceanic state exerts such a strong control over land temperatures (11, 12) as to make the analysis of the small residual difficult. Although a few studies indicated that the fast response of the global or Northern Hemisphere (NH) land temperature to anthropogenic forcings might be notable (1315), an early attempt based on an AGCM found no detectable signal in the spatial pattern of land temperature changes over the second half of the 20th century (16).

Here, we use a series of AGCM simulations to substantiate that the fast response of land temperature to anthropogenic forcings can be detected for NH as a whole in the recent decades. We further show that it is feasible to discern the fast response to aerosols at continental scales, hinting at an emergent constraint on aerosol forcing. We force three Geophysical Fluid Dynamics Laboratory (GFDL) atmospheric models, AM2 (17), AM3 (18), and AM4 (19, 20) with the observed monthly mean sea surface temperatures and sea ice concentrations from 1870 to 2015 (21) but with different forcing combinations (see Materials and Methods). All time-varying forcings [well-mixed greenhouse gases (WMGGs), ozone, aerosols, land use, volcanoes, and solar variations] are held fixed at the 1860 levels in the NO_F experiments. Only aerosols are fixed in NO_AERO, while all forcings, aerosols included, are allowed to vary with time in ALL_F. The difference between NO_AERO and NO_F (the former minus the latter, our convention throughout the paper) is the fast response to nonaerosol forcings, and the difference between ALL_F and NO_AERO is the fast response to aerosol forcing.


Figure 1 compares the simulated NH summer and autumn [June to November (JJASON)] LSAT with observations (22). We focus on these seasons to avoid the large internal atmospheric variability in winter and spring (23). Amid the substantial variations caused by El Niño–Southern Oscillation and major volcanic eruptions, the observed LSAT shows a pronounced multidecadal warming trend since the early 1980s. While capturing the general trend, NO_F does not warm up enough, with a clear departure from the observations after 2000 (Fig. 1A). None of the ensemble members in NO_F can reproduce the observed warming. On average, the NH land warming in the last several decades (defined as the 2001 to 2015 average minus the 1961 to 1980 average) in NO_F is about 77% of the observed value (Table 1). The inclusion of nonaerosol forcings (NO_AERO) allows the models to reproduce the full scale of the observed warming, while the additional contribution from aerosols (i.e., ALL_F versus NO_AERO) appears to have no significant effect (Fig. 1, B and C). This means that about 23% of the observed recent warming can be attributed to the fast response to nonaerosol forcings, but for NH land as a whole, the average fast response to aerosols is negligible in the recent decades. This fast response is presumed to be driven primarily by WMGGs and ozone, rather than land use, volcanoes, or solar variations, whose forcings are relatively small over the past 50 years (24). However, a detailed breakdown will not be needed in the following.

Fig. 1 Time series of NH LSAT anomalies.

Time series of the observed and simulated NH JJASON LSAT anomalies. The reference period is 1961 to 1980. The panels correspond to the three experiments: NO_F (A), NO_AERO (B) and ALL_F (C). The black lines represent the observations. The colored solid lines and shadings are for the ensemble means and ranges, respectively (AM2 in red, AM3 in blue, and AM4 in green). Strong El Niño (E), La Niña (L), and volcanic eruptions (V) are marked. The right panels are blowups of the period 2001–2015 denoted by the rectangles in the left panels. The dashed lines represent the 2001–2015 averages.

Table 1 NH and regional LSAT changes.

Comparison of the observed and simulated recent JJASON LSAT changes (K, 2001–2015 minus 1961–1980) over NH, Europe, and Asia.

View this table:

The fast response to aerosols, which can be estimated by differencing NO_AERO and the observations, becomes evident at continental scales (fig. S1). The responses over Europe and Asia are of opposite signs, and this distinct spatial structure is consistent among the models. The observed and simulated recent LSAT changes (2001–2015 minus 1961–1980) over Europe (30°N to 60°N, 10°W to 50°E) and Asia (5°N to 35°N, 70°E to 120°E) are shown in Fig. 2 (A and B, respectively). As is the case for the entire NH land, the simulated warming in NO_AERO is stronger than in NO_F for both regions, affirming the roles of nonaerosol forcings. NO_AERO, however, significantly underestimates the observed warming by ~25% over Europe but overestimates by ~34% over Asia (Table 1), implying that one may have to invoke the fast response to aerosol forcing to reconcile the model simulations with observations. It is well known that the anthropogenic SO2 emissions (the main source of sulfate aerosol) have declined steadily over Europe since the 1980s owing to air pollution control, while the emissions over Asia (China and India in particular) have substantially increased in the same time period due to economic development (25). It is encouraging that the empirically inferred fast aerosol effects on temperature over Europe and Asia are consistent, at least in sign, with the regional emission trends. This is especially clear in the ALL_F simulations performed with AM3 and AM4, which are appreciably warmer than NO_AERO over Europe and colder over Asia. As a result, the warming in these ALL_F simulations is in good agreement with the observations. In comparison, the AM2-simulated ALL_F does not differ much from its NO_AERO counterpart, suggesting that the aerosol effect in AM2 is much weaker than in the other two models.

Fig. 2 Regional LSAT changes.

Observed and simulated recent JJASON LSAT changes (2001–2015 minus 1961–1980) over Europe (A) and Asia (B). The boxes denote the mean changes (the center lines) and 90% confidence intervals (see Materials and Methods). The observations are in black, AM2 in red, AM3 in blue, and AM4 in green. The black open circles represent the individual ensemble members.

The three models differ substantially in the representation of aerosol physics and chemistry. In AM2, aerosol concentrations are prescribed, and only the direct effects of aerosols on radiation are included, while AM3 and AM4 feature interactive aerosols simulated from prescribed emissions of aerosol precursors and both the direct and indirect effects (cloud albedo and cloud lifetime effects). The global annual mean preindustrial–to–present day (1990) aerosol effective radiative forcing (ERF), computed as in (26), is −0.33 W m−2 in AM2, −1.69 W m−2 in AM3, and − 0.96 W m−2 in AM4 (19), spanning the Intergovernmental Panel on Climate Change Fifth Assessment Report (AR5) range (−0.1 to −1.9 W m−2) (24), while the nonaerosol forcings are largely consistent across the models. The regional JJASON aerosol forcing in AM2 is also much weaker than in the other two models over Europe and Asia (Fig. 3). AM3 and AM4 yield very similar aerosol forcing over Europe but differ by a factor of ~2 over Asia; the direct effects dominate in the former case and the indirect effects in the latter, as the difference in aerosol ERF between AM3 and AM4 is mostly in the indirect effects (19). It is possible that the differences in Asian monsoon (fig. S2) may also contribute to the intermodel spread over Asia.

Fig. 3 Time series of regional ERF.

Time series of the simulated JJASON nonaerosol (the solid lines) and aerosol (the dashed lines) ERF averaged over Europe (A) and Asia (B). AM2 is in red, AM3 in blue, and AM4 in green. The thin and thick lines are for the annual and 5-year running means, respectively.

One can estimate the fast components of the recent LSAT changes by differencing NO_AERO or ALL_F with NO_F (the slow components). As shown in Fig. 4 (A and B), they are positively correlated with the regional ERF in the same period (2001–2015 minus 1961–1980) as more incoming radiation tends to warm the surface. The correlation is weaker over Asia (r2 = 0.74) than over Europe (r2 = 0.92), possibly due to the large natural variability and intermodel differences associated with the Asian monsoon (fig. S2) and the differences in land models. On the other hand, the slopes are very similar between the two regions. Using the observationally (22) based estimates of recent changes in LSAT (2001–2015 minus 1961–1980; 0.55 ± 0.12 K over Europe and −0.05 ± 0.08 K over Asia), one can infer the regional ERF over this time frame as 6.1 ± 1.9 W m−2 over Europe and −1.7 ± 1.9 W m−2 over Asia. The latter implies that anthropogenic aerosols dominate nonaerosol forcings in determining the sign of the net ERF over Asia.

Fig. 4 Emergent constraints for regional and global ERF.

Scatterplots of the fast recent JJASON LSAT changes versus regional ERF (A and B) or global ERF (C and D) over Europe (A and C) and Asia (B and D). AM2 is in red, AM3 in blue, and AM4 in green. The triangles and squares represent NO_AERO and ALL_F (after subtracting the corresponding NO_F), respectively. The best linear fits and prediction errors are represented by the black solid and dashed lines, respectively. The vertical purple lines denote the observational estimates (obtained by subtracting the multimodel average LSAT changes in NO_F from the observations): the mean (the solid lines) and the mean plus and minus one SD (the dashed lines). The gray areas delineate the ranges of the inferred regional or global ERF.

The fast LSAT changes also show strong correlation across these three models with the global annual mean ERF over the entire historical period (2001–2015 minus the 1860 levels; Fig. 4, C and D). The slopes of the best linear fits are of opposite signs over the two regions since the aerosol-induced recent LSAT change is a net warming over Europe but a net cooling over Asia. Nonetheless, they give rise to similar best estimates of the global ERF (1.5 ± 0.7 W m−2 based on the linear fit over Europe and 1.8 ± 0.7 W m−2 based on the linear fit over Asia). This convergence of these observational constraints, which lead to vastly different regional ERF, lends support to the robustness of the inferred global ERF. In the following analysis, we use the total ERF estimate based on LSAT change over Europe because the observational datasets are more robust (see Materials and Methods), and the estimate is less sensitive to the choice of months (see the Supplementary Materials).

One can further estimate the global aerosol forcing at −1.4 ± 0.7 W m−2 by subtracting the nonaerosol forcings (2.9 ± 0.2 W m−2 in the three models) from the inferred total forcing (1.5 ± 0.7 W m−2). This value is appreciably stronger than the best AR5 estimate (−0.9 W m−2) but well within the 90% confidence interval (−0.1 to −1.9 W m−2). It is also within the 68% confidence interval of −0.65 to −1.60 W m−2 provided by (27). The best estimate is at the lowest end of the Coupled Model Intercomparison Project Phase 6 (CMIP6) aerosol ERF range (−0.63 to −1.37 W m−2) (28).

The transient climate response (TCR; defined as the surface temperature change in response to a 1% per year increase of CO2 at the time of doubling), a quantity crucial for near-term climate projection, can be calculated from the historical warming (δT, 0.80 K) (29) and ERF (F) as F2XδT/F, where F2X is the ERF of CO2 doubling (3.8 W m−2) (29). At a historical forcing of 1.5 W m−2 as estimated here, the implied TCR is 2.0 K. This is at the higher end of the AR5 likely range of 1 to 2.5 K (30) but is close to the median TCR of 1.95 K based on CMIP6 models (31).


There have been other attempts to narrow down the large uncertainty of forward model estimates of aerosol forcing by searching for emergent constraints or potentially observable aspects of present or past climates that are correlated across an ensemble of models with the simulated aerosol forcing. Examples are the sensitivity of rain efficiency to aerosols (32), the decadal trend in surface solar radiation (33), and the liquid cloud water response to volcanic eruptions (34) and anthropogenic pollution sources (35). These studies target certain key processes involved in aerosol-cloud interactions but have not been used to assess the total aerosol forcing, with the notable exception of (33), where the 1990 to 2005 solar brightening trend over Europe is used to infer a global mean aerosol forcing of −1.30 ± 0.4 W m−2. Previous studies have estimated aerosol ERF by comparing CMIP5 model simulated and observed temperature change in the 20th century. The best estimates of aerosol ERF in this study is stronger than the estimates in most studies [−0.8 to −1.5 W m−2 in (36), −0.92 W m−2 in (37), and −1.0 W m−2 in (38)], with the exception of (39) (−2.0 W m−2). Note that all these studies are based on AOGCMs, with the limitations discussed in the introduction.

This study takes a completely different approach by providing a framework for constraining the historical aerosol forcing with the observed land surface temperature record using AGCMs. It also empirically shows that the fast response to aerosols is important for determining regional land temperature trends. A major advantage over traditional AOGCM-based detection and attribution studies is that one can avoid the compensating errors in the aerosol forcing and climate sensitivity (40). For example, when coupled to ocean models, AM3 and AM4 have much higher climate sensitivity than AM2. Yet, all three models produce similar LSAT changes in response to oceanic warming and nonaerosol forcings (NO_F and NO_AERO). This is understandable as the atmosphere-land adjustment involves a substantially reduced set of processes as compared to the fully coupled response and offers a solid foundation for reconciling the difference between the observations and NO_AERO simulations by invoking aerosols. By applying AGCM simulations to distinguish the fast and slow responses to anthropogenic forcings and discern the region-specific aerosol effects, we hope to motivate their wider use in understanding regional climate change.

These results are suggestive but limited because of the small number of models used. The uncertainty range of the inferred aerosol forcing might be substantially reduced through an emergent constraint analysis such as that in Fig. 4 if more AGCM simulations driven by different forcing combinations were available, an exercise that can be facilitated by the CMIP6 (41).


Climate data and climate model simulations

We use the observed surface air temperature from the Berkeley Earth Surface Temperatures (BEST), an interpolated dataset based on many preexisting surface air temperature records (22). We also examine other observational datasets, including GISS Surface Temperature Analysis (GISTEMP) (42) and Cowtan and Way (43). The different datasets agree well on LSAT changes over NH and Europe but diverge somewhat over Asia, providing additional motivation for focusing on Europe in our discussion of the emergent constraint on global ERF. The historical simulations are conducted using GFDL AM2, AM3, and AM4, the atmospheric components of three GFDL coupled climate models (CM2.1, CM3, and CM4.0, respectively).

AM2 uses a finite-volume, latitude-longitude dynamical core with a horizontal resolution of ~2° × 2.5° and 24 vertical levels from the surface to ~3 hPa. The shortwave and longwave radiative transfer schemes are as described in (44, 45). Moist convection is parameterized using the relaxed Arakawa-Schubert convection scheme (46). Stratiform clouds are prognosed based on (47), and the treatment of cloud microphysics follows (48, 49). Convective planetary boundary layers are parameterized based on (50), and conventional stability functions dependent on the Richardson number are used for stable boundary layers. Ozone and aerosol concentrations are prescribed in AM2, and only the aerosol direct effects are included.

AM3 uses a finite-volume, cubed-sphere dynamical core and has a horizontal resolution of ~200 km with 48 vertical layers from the surface to ~1 Pa. AM3 implements the parameterizations of deep and shallow cumulus convection as described in (51, 52), and a physically based representation of aerosol-cloud interactions based on a first principles based the parameterization of droplet activation (53). Subgrid variability of vertical velocities in stratiform clouds is parameterized for droplet activation (54). Ozone and aerosols are simulated interactively from emissions with explicit dry and wet removal.

AM4 uses the same finite-volume, cubed-sphere dynamical core as AM3 and has a horizontal resolution of ~100 km with 33 vertical layers from the surface to ~1 hPa. It includes a double-plume scheme representing deep and shallow convection. The stratiform cloud and boundary layer schemes in AM4 are essentially the same as in AM3. AM4 includes a simplified chemistry scheme, which makes it possible to simulate aerosols interactively, but ozone is prescribed. The cloud microphysics and aerosol schemes in AM4 are similar to those in AM3 but with substantial changes in droplet activation and wet removal of aerosols by convection and snow/ice.

The treatment of aerosol forcing in AM3 and AM4 differs substantially from that in AM2; AM3 and AM4 include both the direct and indirect effects and treat sulfate and hydrophilic black carbon as internally mixed for the direct effects. The formulations and performances of AM2, AM3, and AM4 have been documented in previous papers (1720).

As in Atmospheric Model Intercomparison Project protocols, all simulations in this study are forced with the observed monthly sea surface temperatures and sea ice concentrations (21). ALL_F is driven by the historical radiative forcing agents used for CMIP5 (55) for 1870–2005 and by the Representative Concentration Pathway of 4.5 for 2006–2015. The forcing agents include WMGGs, ozone, aerosols, land use, volcanoes, and solar variations. They are all fixed at the 1860 levels (preindustrial conditions) in NO_F, and only aerosols are fixed in NO_AERO. NO_F and NO_AERO are run from 1870 to 2015 as well. It has been shown that using CMIP6 emissions yields a similar global aerosol ERF (20). For aerosols, AM2 uses prescribed concentrations (either fixed or time varying), while AM3 and AM4 use prescribed emissions and compute their atmospheric concentrations online. Each simulation consists of five ensemble members with different initial conditions to account for internal atmospheric variability.

Statistical analysis

Here, the recent LSAT change is defined as the 2001–2015 average minus the 1961–1980 average. The Student’s t test is applied to assess the statistical significance and confidence intervals of the LSAT changes in the observational record and model ensemble means in Fig. 2 and figs. S1A and S2. The t statistic is calculated ast=X¯2X¯1Sp1n1+1n2where X is the mean surface air temperature and n is the number of years. Subscripts 1 and 2 denote the time periods 1961–1980 and 2001–2015, respectively. Sp is the pooled estimate of the SDSp=(n11)S12+(n21)S22n1+n21where S2 is the unbiased estimate of the variance of surface air temperature. The degree of freedom of the t test is n1 + n2 – 2. Least-squares linear regressions and prediction errors in Fig. 4 and fig. S3 are calculated using the ordinary least-squares regression with Python statsmodel package.

In fig. S1, ensemble members are first interpolated linearly to the observational grid and then averaged for comparison with the observations. A linear contrast test (56) is conducted to assess the statistical significance of the difference between observations and model ensemble means in Table 1 and fig. S1. In the contrast test, the F statistic is calculated asF=X¯4X¯3X¯2+X¯1MSW1n1+1n2+1n3+1n4where X1 and X2 are the 1961 to 1980 and 2001 to 2015 average temperature from observations, respectively, and X3 and X4 are the 1961 to 1980 and 2001 to 2015 average temperature from simulations, respectively. Mean square within (MSW) is the mean square error within the four groups from the one-way analysis of variance (ANOVA), and n1, n2, n3, and n4 are the number of years in each group (n1 = 20, n2 = 15, n3 = 20, and n4 = 15). The calculated F statistic is then compared with the F distribution with (1, n1 + n2 + n3 + n4 – 4) degrees of freedom.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank T. Knutson and D. Paynter for helpful comments on an earlier draft. Funding: Z.S. is funded during the early stage of the project by the National Oceanic and Atmospheric Administration and U.S. Department of Commerce under award NA14OAR4320106. The statements, findings, conclusions, and recommendations are those of the authors and do not necessarily reflect the views of the National Oceanic and Atmospheric Administration or the U.S. Department of Commerce. Author contributions: All authors collaboratively designed the study and wrote the manuscript. Z.S. carried out the model simulations and conducted the analysis. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. BEST was downloaded from the Berkeley Earth website ( Model simulations are archived at GFDL and are available from the authors upon request. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article