Articles | Volume 23, issue 20
Research article
 | Highlight paper
18 Oct 2023
Research article | Highlight paper |  | 18 Oct 2023

Global observations of aerosol indirect effects from marine liquid clouds

Casey J. Wall, Trude Storelvmo, and Anna Possner

Interactions between aerosols and liquid clouds are one of the largest sources of uncertainty in the historical radiative forcing of climate. One widely shared goal to reduce this uncertainty is to decompose radiative anomalies arising from aerosol–cloud interactions into components associated with changes in cloud-droplet number concentration (Twomey effect), liquid-water-path adjustments, and cloud-fraction adjustments. However, there has not been a quantitative foundation for simultaneously estimating these components with global satellite observations. Here we present a method for assessing shortwave radiative flux anomalies from the Twomey effect and cloud adjustments over ocean between 55 S and 55 N. We find that larger aerosol concentrations are associated with widespread cloud brightening from the Twomey effect, a positive radiative adjustment from decreasing liquid water path in subtropical stratocumulus regions, and a negative radiative adjustment from increasing cloud fraction in the subtropics and midlatitudes. The Twomey effect and total cloud adjustment have contributed 0.77 ± 0.25 and 1.02 ± 0.43 W m−2, respectively, to the effective radiative forcing since 1850 over the domain (95 % confidence). Our findings reduce uncertainty in these components of aerosol forcing and suggest that cloud adjustments make a larger contribution to the forcing than is commonly believed.

1 Introduction

Changes in aerosol concentrations over the industrial era have modified clouds and perturbed the global radiation balance at the top of the atmosphere (Raghuraman et al., 2021; Kramer et al., 2021). The radiative flux perturbation resulting from these cloud changes, known as the effective radiative forcing from aerosol–cloud interactions (ERFaci), is estimated to be 0.84 ± 0.61 W m−2 between 1750 and 2019 (90 % confidence interval, CI, from Forster et al., 2021). ERFaci is much more uncertain than the positive radiative forcing from carbon dioxide changes (+2.16 ± 0.26 W m−2), meaning that ERFaci offsets a potentially large but highly uncertain portion of historical greenhouse-gas forcing. Reducing this uncertainty would improve assessments of climate sensitivity and committed future warming (Matthews and Zickfeld, 2012; Mauritsen and Pincus, 2017; Sherwood et al., 2020; Watson-Parris and Smith, 2022).

An extension of these forcing estimates is to characterize how changes in different cloud properties contribute to ERFaci, thereby providing insight into the relative importance of different processes. For instance, as cloud condensation nuclei (CCN) become more abundant, liquid clouds typically form smaller but more numerous droplets. The change in cloud droplet effective radius (re) and number concentration (Nd) directly increases cloud optical thickness – a mechanism known as the Twomey effect (Twomey, 1977). The reduction of cloud droplet size can also enhance evaporation or reduce precipitation, causing adjustments in cloud thickness, lifetime, or morphology (Albrecht, 1989; Pincus and Baker, 1994; Rosenfeld et al., 2006; Bretherton et al., 2007). Separating the radiative impacts of the Twomey effect and cloud adjustments is thus an important step towards understanding the causes of ERFaci.

Recent community assessments find that the components of ERFaci all have considerable uncertainty. The Sixth Assessment Report of the Intergovernmental Panel on Climate Change (IPCC) estimates that the Twomey effect is 0.7 ± 0.5 W m−2, the adjustment of liquid water path (LWP) is +0.2 ± 0.2 W m−2, and the adjustment of liquid-cloud fraction is 0.5 ± 0.4 W m−2 (90 % CIs for forcing between 1750 and 2014) (Forster et al., 2021). Another assessment by the World Climate Research Programme reports even larger uncertainties (Bellouin et al., 2020). In particular, they find that the cloud-fraction adjustment is especially uncertain: it could be negligibly small or large enough to offset most of the historical carbon-dioxide forcing.

Constraints from satellite observations offer a path toward reducing this uncertainty, but in practice it has been difficult to isolate relationships between aerosols and radiative anomalies caused by changes in individual cloud properties (Feingold et al., 2022). Furthermore, previous global observational studies that attempt to quantify these relationships used separate methods for estimating the Twomey effect, LWP adjustment, and cloud-fraction adjustment, so their estimates may suffer from limitations that differ from one ERFaci component to another (Forster et al., 2021). This complicates efforts to rigorously compare the Twomey effect and cloud adjustments. Here we address these challenges by adapting techniques from the cloud-feedback literature. We develop a cloud radiative kernel that separates the radiative anomalies caused by changes in re, LWP, and liquid-cloud amount, and we relate each of these radiative anomalies to local aerosol concentrations and Nd. This facilitates an assessment of the Twomey effect and liquid-cloud adjustments over the global ocean.

2 Data and methods

2.1 Satellite data, reanalysis, and climate model output

We analyze monthly gridded satellite observations from 2003 through 2020 obtained from the Moderate Imaging Spectroradiometer (MODIS) MCD06COSP dataset version 6.2.0 (Pincus et al., 2023). This dataset combines observations from MODIS instruments on the Aqua and Terra satellites. Our primary unit of analysis is a joint histogram of pixel counts for liquid-topped clouds partitioned by re and LWP. Histogram counts are normalized by the number of all valid pixels in the grid box and then multiplied by 100 to convert the units to cloud fraction (Fig. 1a). These histograms represent the fractional occurrence of liquid-topped clouds that are exposed to space; hence they do not include cases where liquid cloud is obscured by overlying ice. LWP is defined by the vertical integral of cloud liquid water mass per unit area, and re is defined by the ratio of the third and second moments of the cloud-droplet radius distribution. These two variables are estimated with the MODIS 3.7 µm retrieval algorithm (Platnick et al., 2017).

Figure 1MODIS joint histogram and SW cloud radiative kernel averaged over the latitude, longitude, and time dimensions. Averages are computed over ocean between 55 S and 55 N. (a) Joint histogram of liquid-cloud fraction (C) partitioned by liquid water path (LWP) and cloud-droplet effective radius (re). (b) SW cloud radiative kernel. The kernel represents R/C, where R is the SW radiative effect of liquid clouds at the top of the atmosphere.


We also use daily gridded Nd estimates from MODIS (Gryspeerdt et al., 2022a) and monthly gridded radiative flux retrievals from the Clouds and the Earth's Radiant Energy System (CERES) Energy Balanced and Filled ed. 4.1 and FluxByCldTyp ed. 4.1 datasets (Loeb et al., 2018; Sun et al., 2022). The MODIS Nd estimates can be biased when re is sufficiently small, when the cloud visible optical thickness is sufficiently small, or when three-dimensional radiative transfer effects contribute to the measured radiances (Grosvenor et al., 2018). To avoid these problematic cases, Nd is estimated for single-layer liquid clouds that satisfy several conditions: (i) re is larger than 4 µm; (ii) cloud optical thickness is larger than 4; (iii) cloud fraction at 5 km resolution is larger than 0.9; (iv) the solar zenith angle is less than 65; (v) the satellite viewing zenith angle is less than 55; (vi) the sub-pixel heterogeneity index defined by Zhang and Platnick (2011) is less than 0.3; (vii) the MODIS estimate of re is largest from the 3.7 µm retrieval algorithm, followed by the 2.1 µm retrieval algorithm and then the 1.6 µm retrieval algorithm; and (viii) cloud optical thickness is in the top 10 % of values in 100 km× 100 km regions (“Z18 sampling” in Gryspeerdt et al., 2022a). The final condition preferentially selects the convective cores in cloudy scenes (Zhu et al., 2018). Most cloud droplets in shallow convective clouds form near the cloud base, so Nd in convective cores depends on the CCN concentration in air that enters the cloud from below (Rosenfeld et al., 2019). Thus, Nd in convective cores typically does not represent Nd in the entire cloud, but it serves as an indicator of CCN concentration near cloud base. We denote the estimated cloud-droplet number concentration as Ñd to distinguish it from the cloud-droplet number concentration in the entire cloud. For consistency with the MODIS cloud histograms, we use Ñd estimates from the MODIS 3.7 µm retrieval algorithm, and we combine data from the Terra and Aqua satellites. Daily values of lnÑd are averaged across the satellite platforms and over 1-month intervals, weighted by the number of pixels with a valid Ñd retrieval.

Monthly meteorological fields and the dry-mass concentration of sulfate aerosol at 910 hPas, are obtained from MERRA-2 reanalysis (Gelaro et al., 2017; Randles et al., 2017). We consider sulfate aerosol because it dominates the anthropogenic influence on CCN (Charlson et al., 1992; Stevens, 2015), and we select data from 910 hPa rather than the surface because the 910 hPa level is a better indicator of aerosol concentration near cloud base (Painemal et al., 2017). The sulfate data are determined from bias-corrected observations of total aerosol optical depth from cloud-free pixels and simulations from a global model that represents the sources, sinks, and chemistry of sulfate and its precursor gases. The data assimilation accounts for aerosol swelling in humid environments and filters out pixels near clouds that are affected by retrieval bias (Randles et al., 2017). The main limitation of the sulfate data is that the total aerosol optical depth is constrained by observations, but aerosol species distributions and vertical profiles are not. These data provide an additional indicator of cloud-base CCN concentration that is independent of the MODIS estimates of Ñd.

Finally, we use output from historical simulations of 20 global climate models (GCMs) from the Coupled Model Intercomparison Project Phase 6 (CMIP6). The simulations are run from 1850 through 2014 with realistic emissions of greenhouse gases, aerosols, and aerosol precursor gases. Sulfate mass concentration from the model output is converted to pressure coordinates and linearly interpolated to 910 hPa. The models are listed in Table S1 in the Supplement.

2.2 Quantifying aerosol indirect effects

We first relate variability of aerosols to the radiative effects of liquid clouds – relationships we call “aerosol indirect effects.” Our analysis begins with the MODIS joint histograms of liquid-cloud fraction, C, partitioned by re and LWP (Fig. 1a). For a given latitude, longitude, and time, the shortwave (SW) radiative flux anomaly at the top of the atmosphere that is induced by liquid clouds, R, is estimated according to


where r and l represent the re and LWP dimensions of the histogram, and primes denote monthly anomalies relative to the local climatological seasonal cycle. The cloud radiative kernel, R/Crl, represents the SW flux anomaly that would occur if the cloud fraction Crl were to increase by 1 % with all non-cloud factors fixed (Fig. 1b). The kernel is computed with the RRTMG radiative transfer model (Clough et al., 2005) following a method similar to that of Zelinka et al. (2012). We then adapt a method of cloud-feedback analysis developed by Zelinka et al. (2013) to decompose R into contributions from different cloud properties:


where Rre, RLWP, and RCF are the SW flux anomalies caused by re, LWP, and cloud-amount anomalies, respectively – each computed with the other properties held fixed. We note that Rre is the radiative anomaly that is caused by variations in re with fixed LWP, so it is equivalent to the radiative anomaly that is caused by variations in Nd with fixed liquid-water content. Rres is the residual of the decomposition. The radiative kernel and MODIS joint histograms reproduce monthly observations of R across the global ocean with a bias of about +4.6 % (Appendix A). The methods for computing the kernel and decomposing R are described in Appendix A and B.

To test robustness of the results, we make one set of R estimates in which only fully cloud-covered pixels are included in the histograms and a second set of estimates in which fully and partly cloud-covered pixels are both included. We refer to these cases as MODISCLD and MODISCLD+PCL, respectively. The filter of the MODISCLD case avoids retrieval biases that affect partly cloudy pixels, but it may introduce a sampling bias by excluding some cloud elements. The opposite is true for MODISCLD+PCL. Both cases are presented to explore trade-offs between the accuracy and completeness of the satellite cloud data.

We relate R to sulfate and cloud-droplet concentrations using cloud-controlling factor analysis (Scott et al., 2020; Myers et al., 2021). Our analysis closely follows the method of Wall et al. (2022; hereafter W22) except that we generalize their results by applying cloud-controlling factor analysis to R and each of its components. The cloud-controlling factor method approximates R as a linear combination of seven local cloud-controlling factors xi:


The first six xi terms include sea surface temperature, estimated inversion strength at the top of the planetary boundary layer (Wood and Bretherton, 2006), low-level advection across a surface-temperature gradient, surface wind speed, relative humidity at 700 hPa, and vertical wind at 700 hPa. Collectively these terms represent all of the standard large-scale meteorological controls on liquid clouds in the marine boundary layer that have been proposed in the literature (Scott et al., 2020). The final xi term can be either ln s or lnÑd. All meteorological terms and ln s are calculated with MERRA-2 data and linearly interpolated to the native 1× 1 grid of MODIS. We then select ocean-covered grid boxes, remove the climatological seasonal cycle and least-squares linear trend from all variables in each grid box, and average the anomalies over a 5× 5 grid that spans 55 S to 55 N. For each ocean grid box, R is regressed against anomalies of the seven cloud-controlling factors using ordinary least-squares multilinear regression. Separate regressions are performed with ln s and lnÑd as the final predictor. Thus, the regression coefficients R/lns and R/lnÑd represent the relationship between R and local anomalies of ln s or lnÑd with all meteorological predictors held constant. The goodness of fit of the regression model is determined by computing the fraction of R variance that it explains in each grid box and then spatially averaging the results over the domain. On average, the regression method explains 46 % of the R variance when ln s is the final predictor and 49 % of the variance when lnÑd is the final predictor. We also apply the method to Rre, RLWP, and RCF to estimate the Twomey effect, LWP adjustment, and cloud-fraction adjustment.

Our analysis differs from existing global estimates of aerosol indirect effects in several ways. First, estimates that control for fewer meteorological factors are susceptible to bias from correlations between meteorology and aerosols (Mauger and Norris, 2007; Gryspeerdt et al., 2016; Andersen et al., 2017). Our method minimizes this bias by controlling for all of the standard large-scale meteorological drivers of liquid boundary-layer clouds that have been proposed in the literature. Second, estimates derived from daily or monthly grid-box averages of re and LWP suffer from aggregation bias because these properties are nonlinearly related to cloud albedo (Feingold et al., 2022; Goren et al., 2023). Our method avoids this bias by inferring cloud radiative effects from joint histograms of re and LWP rather than grid-box averages. Third, studies of ship tracks, industrial plumes, or volcanic eruptions offer some of the most captivating evidence of aerosol indirect effects, but the estimates they provide may not be representative of the global scale (Possner et al., 2018; Toll et al., 2019; Glassmeier et al., 2021). Our method avoids this potential sampling bias by estimating aerosol indirect effects across the global ocean. Fourth, no global observational study has simultaneously estimated the Twomey effect, LWP adjustment, and cloud-fraction adjustment, so comparisons of these components have been complicated by the fact that each one is estimated using different data, methods, assumptions, and uncertainty quantification (Forster et al., 2021). Our method avoids this complication by estimating all components with a single, self-consistent framework. Our estimates of aerosol indirect effects could still be affected by satellite retrieval biases, but they improve upon existing estimates in these four ways (Painemal and Zuidema, 2011; Ma et al., 2018).

3 Global analysis of aerosol indirect effects

We next present estimates of the Twomey effect and cloud adjustments across the global ocean. The Twomey effect can occur whenever the CCN concentration is small enough that it limits the number of droplets that form in cloud updrafts. This condition is usually satisfied over ocean, so the Twomey effect is expected to be ubiquitous in oceanic clouds (Rosenfeld et al., 2014). Indeed, we find that increasing sulfate concentration is associated with significant cloud brightening from Rre changes across most of the global ocean (Fig. 2a). Cloud brightening is also observed in response to increasing Ñd (Fig. 2b), but we caution against overinterpreting statistical significance of this relationship because Rre and Ñd are both inferred from the MODIS re retrievals (Ñd is inferred using 10 % of the data). Nevertheless, these results show that most marine liquid clouds exhibit the Twomey effect.

Figure 2Aerosol indirect effects estimated with two indicators of cloud-base CCN concentration: sulfate aerosol mass concentration at 910 hPa (s) and cloud-droplet number concentration from pixels with the largest 10 % cloud optical thickness (Ñd). Linear regression coefficients are shown for (a) Rre/lns, (b) Rre/lnÑd, (c) RLWP/lns, (d) RLWP/lnÑd, (e) RCF/lns, and (f) RCF/lnÑd, where Rre, RLWP, and RCF are the top-of-atmosphere SW flux perturbations caused by re anomalies, LWP anomalies, and cloud-fraction anomalies, respectively. Panels (a, b) represent the Twomey effect, panels (c, d) represent the LWP adjustment, and panels (e, f) represent the cloud-fraction adjustment. Stippling indicates regression coefficients that are significantly different from zero with the false discovery rate limited to 0.1 (Wilks, 2016). Cloud radiative effects are computed with only fully cloud-covered pixels included in the cloud histograms (MODISCLD). Note that the contour values in (a), (c), and (e) are proportional to those in (b), (d), and (f).

In contrast, previous work has shown that cloud adjustments can differ from one cloud regime to another (Zhang et al., 2022). As cloud droplets become smaller, they sediment more slowly out of the cloud-top entrainment zone, and they evaporate more quickly when exposed to entrained air. This enhances evaporation and reduces LWP in non-precipitating clouds (Bretherton et al., 2007; Small et al., 2009). Clouds with smaller droplets also form precipitation more slowly through collision and coalescence. This may cause other changes in cloud properties, including deeper cumulus clouds, longer cloud lifetimes, larger stratiform areas detrained from precipitating cloud elements, or changes in mesoscale cellular structure (Albrecht, 1989; Pincus and Baker, 1994; Rosenfeld et al., 2006; Seifert et al., 2015; Possner et al., 2018; Dagan et al., 2017; Goren et al., 2022). The cloud adjustments can, in turn, affect CCN concentration by changing precipitation scavenging or sulfate formation in cloud droplets (Wood et al., 2012; Kang et al., 2022; Andreae and Rosenfeld, 2008). Some of these mechanisms depend on the meteorological conditions, so they may vary regionally (Chen et al., 2014; Possner et al., 2020; Zhou et al., 2021; Zhang and Feingold, 2023).

The estimated cloud adjustments exhibit regional variations that are consistent with some of these proposed mechanisms. The radiative adjustment from LWP changes is positive in much of the subtropics, and it maximizes in areas of semi-permanent stratocumulus clouds and directly downwind (Fig. 2c and d). This spatial pattern suggests that enhanced evaporation in non-precipitating stratocumulus may contribute to the LWP adjustment (Bender et al., 2019). Weak or insignificant LWP adjustments are found across much of the midlatitude oceans, despite the fact that precipitating clouds occur less frequently in these regions as sulfate aerosols become more abundant (W22). This weak overall LWP adjustment might be partly explained by the fact that different cloud regimes in extratropical cyclones exhibit adjustments that may counteract one another (Naud et al., 2017; McCoy et al., 2018). Furthermore, the radiative adjustment from cloud-fraction changes is negative in most of the subtropics and midlatitudes, suggesting that adjustments in these regions may involve changes in cloud lifetime, size, or morphology as well (Fig. 2e and f). Subtropical stratocumulus regions exhibit offsetting adjustments from LWP changes and cloud-fraction changes. In these cases, as CCN become more abundant, the overall liquid-cloud fraction increases, but the increase is disproportionately large in cloud elements with below-average LWP. This combination is consistent with larger stratiform areas detrained from precipitating clouds or a shift from open to closed mesoscale cellular convection (Possner et al., 2018; Rosenfeld et al., 2006). Thus, aerosol-driven changes in evaporation and precipitation may both contribute to cloud adjustments. The spatial patterns of the adjustments predicted with ln s resemble those predicted with lnÑd, suggesting that the estimated adjustments are robust.

We determine the relative importance of the Twomey effect and cloud adjustments at the global scale by averaging the regression coefficients over ocean. These averages can be interpreted as the cloud radiative anomalies that would occur if sulfate concentration and Ñd were increased by a factor of 2.7 at every location. Uncertainty quantification for these estimates is described in Appendix C. For the regressions against ln s and lnÑd, we find that the Twomey effect and cloud-fraction adjustment both significantly increase SW reflection to space (Fig. 3). The cloud-fraction adjustment enhances SW reflection by between 43 % and 250 % compared to the Twomey effect alone (95 % CI), so it makes a substantial contribution to the overall aerosol indirect effect. The large relative magnitude of the cloud-fraction adjustment is consistent with the observed cloud response during a volcanic eruption in Holuhraun, Iceland, and the uncertainty range for this adjustment overlaps with other observational estimates (Chen et al., 2022; Gryspeerdt et al., 2020). Furthermore, the LWP adjustment reduces SW reflection to space, offsetting between 6 % and 87 % of the Twomey effect. This uncertainty range is comparable to the range of estimates from Diamond et al. (2020) and Gryspeerdt et al. (2019). However, it differs from a ship-track analysis by Manshausen et al. (2022), who find that the LWP adjustment increases SW reflection rather than reduces it. This apparent discrepancy may be a consequence of the fact that localized, short-term aerosol perturbations such as shipping emissions can give rise to different cloud adjustments than sustained aerosol perturbations over larger spatial and temporal scales (Glassmeier et al. 2021). Although the LWP adjustment and cloud-fraction adjustment counteract one another, the total cloud adjustment is still significantly negative, and it is comparable to the Twomey effect. Furthermore, the estimated Twomey effect and total cloud adjustment are similar for the MODISCLD and MODISCLD+PCL cases, indicating that they do not change substantially when partly cloudy pixels are filtered in different ways. The Twomey effect and total cloud adjustment are also qualitatively consistent when different Nd datasets are used, and they are an order of magnitude larger than Rres/lns and Rres/lnÑd (Figs. S1 and S2 in the Supplement). These results robustly show that the Twomey effect and total cloud adjustment cause comparable changes in top-of-atmosphere SW flux.

Figure 3Spatial averages of the regression coefficients that represent aerosol indirect effects. Averages are computed over ocean between 55 S and 55 N. (a) Aerosol indirect effects estimated with s as the CCN indicator. RLWP+CF/lns represents the total cloud adjustment (RLWP+CFRLWP+RCF). The MODISCLD case is computed with only fully cloud-covered pixels included in the cloud histograms, and the MODISCLD+PCL case is computed with both fully and partly cloud-covered pixels included. (b) Similar to (a) except that Ñd is the CCN indicator. Squares show mean values, and vertical lines show 95 % CIs.


A limitation of these results is that the MODISCLD and MODISCLD+PCL cases have offsetting differences in the estimated LWP and cloud-fraction adjustments (Fig. 3). This means that the estimates of the individual LWP and cloud-fraction adjustments depend on filtering of partly cloudy pixels, but the estimate of the total cloud adjustment does not. One implication is that aerosol variations must be associated with changes in the relative amounts of partly and fully cloud-covered pixels. Partly and fully cloudy pixels reside on cloud edges and interiors, respectively, so a change in the relative amounts of these pixels implies a change in the cloud perimeter-to-area ratio (W22). This suggests that the global cloud adjustment may involve changes in the horizontal size or morphology of the clouds. A case study demonstrating this concept is presented in Appendix D. A second implication is that whenever aerosol perturbations change the cloud size or morphology, the conventional practice of estimating the individual LWP and cloud-fraction adjustments will inevitably lead to results that depend on the classification and filtering of partly cloudy pixels. Thus, instrument sensitivity, horizontal resolution, and subjective pixel-classification thresholds can all affect the results. In contrast, we find that more robust results can be obtained by estimating the total cloud adjustment. Our analysis provides the first direct assessment of this quantity from observations, reanalysis, and radiative transfer modeling.

4 Implications for historical aerosol forcing

We next combine estimates of aerosol indirect effects and historical sulfate changes to infer ERFaci. The forcing estimates use sulfate rather than Ñd because sulfate concentration is a widely available output variable from GCMs but Ñd is not. Assuming that sulfate dominates the anthropogenic influence on CCN (Charlson et al., 1992; Stevens, 2015), we estimate SW ERFaci from liquid-topped clouds according to


where Δln s is the change in sulfate concentration between preindustrial (1850–1859) and present-day (2005–2014) conditions simulated by GCMs that participated in CMIP6. This method of estimating ERFaci has been validated with volcanic eruptions and other known variations of regional sulfur-dioxide emissions (W22).

We note that this method differs from that of a similar study by W22 in three ways. First, we estimate ERFaci from all liquid-topped clouds, while W22 estimate ERFaci from low-level clouds, defined as clouds with tops between the surface and 680 hPa. We applied their method to estimate SW ERFaci from all liquid-topped clouds, and we found that the result is about 26 % larger in magnitude than the estimate of SW ERFaci from low-level clouds. Second, we estimate SW ERFaci, while W22 estimate net ERFaci. Thus, their estimate includes an additional ERFaci component from changes in longwave radiation, which offsets about 14 % of the SW component (Appendix B). Third, we estimate ERFaci with MODIS observations and radiative kernels, while W22 estimate ERFaci with CERES observations for their main result. Our kernel-based estimates of R are about 4.6 % larger in magnitude than the corresponding CERES observations (Appendix A). These three factors cause our ERFaci estimates to have a larger magnitude than those reported by W22.

Our method for estimating ERFaci can be applied to each component of the aerosol indirect effect to quantify the associated historical effective radiative forcing. The Twomey effect contributes a negative instantaneous radiative forcing (IRFaci) that peaks in subtropical stratocumulus regions and the midlatitude oceans of the Northern Hemisphere (Fig. 4a). Forcing is relatively large in stratocumulus regions because these clouds exhibit a strong radiative response to CCN perturbations (Fig. 2), and forcing is relatively large in the Northern Hemisphere midlatitudes because these regions are close to anthropogenic aerosol sources. The geographic pattern and magnitude of the IRFaci generally agree with the estimates of McCoy et al. (2017), Kinne (2019), and Jia et al. (2021). Furthermore, the IRFaci is comparable to the effective radiative forcing from the combined effect of LWP adjustments (ALWP) and cloud-fraction adjustments (ACF) (Fig. 4b). Thus, the magnitude of the overall ERFaci peaks in the subtropical stratocumulus regions and the Northern Hemisphere midlatitudes as well (Fig. 4c).

Figure 4Components of historical ERFaci from liquid clouds, including (a) the Twomey effect (IRFaci), (b) the total cloud adjustment (ALWP+ACF), and (c) the overall ERFaci. The estimates represent SW forcing, and they are computed with only fully cloud-covered pixels included in the cloud histograms (MODISCLD).

We average the forcing components over ocean between 55 S and 55 N to quantify their large-scale climate impacts. Confidence intervals are computed accounting for regression-slope uncertainty and inter-model spread in the estimates of Δln s (Appendix C). To frame our results in the context of the existing literature, we compare our estimates with forcing calculations from an assessment of the World Climate Research Programme (WCRP) (Bellouin et al., 2020) and forcing estimates from 14 GCMs that participated in the Coupled Model Intercomparison Project Phase 5 (CMIP5) and AeroCom experiments computed by Gryspeerdt et al. (2020) (Table S1). We repeat the original WCRP analysis except that we restrict the calculation to SW forcing over ocean between 55 S and 55 N, as described in the Supplement. The GCM forcing estimates are averaged over ocean between 55 S and 55 N as well. Although all forcing estimates are computed over the same spatial domain, their time periods differ slightly from one another: The present-day reference years are 2005–2014 for our estimates, 2005–2015 for the WCRP estimates, and 2000 for the CMIP5 and AeroCom estimates, and the preindustrial reference years are 1850–1859 for our estimates, 1850 for the WCRP and CMIP5 estimates, and 1860 for the AeroCom estimates (Bellouin et al., 2020; Zelinka et al., 2014; Ghan et al., 2016). These differences in the reference periods could cause differences in ERFaci of 0.1 W m−2 or less (IPCC, 2021).

Figure 5Components of SW ERFaci for liquid clouds averaged over ocean between 55 S and 55 N. The “MODISCLD” and “MODISCLD+PCL” cases show estimates from this study. The “WCRP” case is computed with the method of Bellouin et al. (2020). The “GCMs” case shows values from 14 CMIP5 and AeroCom models computed by Gryspeerdt et al. (2020). The WCRP and GCM estimates are modified from their original published values so that they represent averages over ocean between 55 S and 55 N. Squares show median values, vertical lines show 95 % CIs, and dots in the GCMs case show individual models.


Averages of IRFaci, cloud adjustments, and the overall ERFaci are compared in Fig. 5. We find that ACF is significantly negative, and ALWP is more likely than not to be positive, but the magnitudes of these estimates depend on filtering choices for partly cloudy pixels. The three remaining components are insensitive to partly cloudy pixel filtering: the IRFaci is 0.77 ± 0.25 W m−2, the total cloud adjustment is 1.02 ± 0.43 W m−2, and the overall ERFaci is 1.86 ± 0.62 W m−2 (95 % CIs from MODISCLD). These results lie inside the ranges of the corresponding WCRP and GCM estimates, and the median IRFaci agrees very well with the WCRP and GCM values. However, our results reduce uncertainty of each component by at least 62 % relative to the confidence intervals of the WCRP and at least 23 % relative to the range of GCMs. Furthermore, the WCRP and GCM estimates do not rule out the possibility that the total cloud adjustment is positive or an order of magnitude smaller than the IRFaci. According to our analysis, however, such a small or positive adjustment is implausible. Thus, our analysis reduces uncertainty in the historical IRFaci and total cloud adjustment, and it clarifies the relative importance of these components.

The ocean-average SW ERFaci can be scaled to establish an upper bound for the global-average net ERFaci. Let ERFaci,net,g be the global-average net ERFaci, and let ERFaci,net,d be the domain-average net ERFaci, where the domain includes ocean areas between 55 S and 55 N. Assuming that the average net ERFaci is negative outside the domain (Diamond et al., 2020), it follows that

(1) ERF aci , net , g < f ERF aci , net , d ,

where f=0.56 is the fraction of global surface area covered by the domain. ERFaci,net,d can be expressed as

(2) ERF aci , net , d = ERF aci , sw , d ( 1 - β ) ,

where ERFaci,sw,d is the domain-average SW ERFaci, and β is the fraction of SW ERFaci that is offset by longwave ERFaci. Our radiative kernel cannot accurately assess longwave cloud radiative effects, so β is estimated by applying cloud-controlling factor analysis to CERES satellite data instead (Appendix B). The resulting uncertainty range is β= 0.14 ± 0.06 (90 % CI). We evaluate Eq. (2) with the bounds of the 90 % CIs of β and ERFaci,sw,d from the MODISCLD and MODISCLD+PCL cases and then select the least-negative value of ERFaci,net,d to evaluate Eq. (1). The result constitutes an upper bound for ERFaci,net,g.

The above reasoning implies a 95 % probability that the global net ERFaci from liquid clouds is more negative than 0.56 W m−2 (relative to 1850–1859). Equivalent upper bounds from the published literature include 0.07 W m−2 from the WCRP assessment and 0.3 W m−2 from the observation-based estimate of the IPCC Sixth Assessment Report (relative to 1850 and 1750, respectively) (Bellouin et al., 2020; Forster et al., 2021). Our analysis thus supports a more stringent upper bound on global ERFaci. This constraint is similar to another estimate from cloud-controlling factor analysis presented by W22, but our estimate invokes weaker assumptions because it does not extrapolate forcing to areas outside the domain. Our upper-bound estimate also complements evidence from global energy-balance arguments, which constrains the lower bound of ERFaci (Stevens, 2015; Smith et al., 2021).

5 Conclusion

We analyze MODIS satellite data and adapt techniques from the cloud-feedback literature to quantify aerosol indirect effects from liquid-topped clouds. Our method avoids aggregation and sampling biases that may affect some previous studies, and it controls for all of the standard large-scale meteorological drivers of liquid boundary-layer clouds that have been proposed in the literature, thereby minimizing biases from confounding meteorological factors (Possner et al., 2018; Glassmeier et al., 2021; Feingold et al., 2022). Furthermore, the Twomey effect, LWP adjustment, and cloud-fraction adjustment are simultaneously quantified with a single, self-consistent framework. This guarantees that all of the components, and their uncertainties, are quantified in a consistent way. Although it is important to continue characterizing satellite retrieval biases, to include new meteorological cloud-controlling factors as they are discovered, to investigate nonlinear and nonlocal relationships between clouds and their controlling factors (Lewis et al., 2023), and to quantify additional ERFaci components from ice-containing clouds, our method overcomes several limitations that affect previous observational estimates of aerosol indirect effects.

We apply our method to constrain aerosol indirect effects across the global ocean. We find that increasing CCN concentration is associated with widespread cloud brightening from the Twomey effect, a positive radiative adjustment from decreasing LWP in subtropical stratocumulus regions, and a negative radiative adjustment from increasing cloud fraction in the subtropics and midlatitudes. The estimated aerosol indirect effects are combined with historical sulfate changes simulated by CMIP6 models to quantify the associated SW ERFaci. The Twomey effect and total cloud adjustment are estimated to contribute 0.77 ± 0.25 and 1.02 ± 0.43 W m−2, respectively, to the SW ERFaci averaged over ocean between 55 S and 55 N (95 % CIs). Our findings reduce uncertainty in these components of aerosol forcing and suggest that liquid-cloud adjustments make a larger contribution to the forcing than is commonly believed.

Appendix A: Cloud radiative kernel

We compute a SW cloud radiative kernel for the MODIS re–LWP joint histogram to quantify the effect of liquid-cloud anomalies on the top-of-atmosphere SW flux. The radiative kernel is similar to that of Zelinka et al. (2012) with two exceptions. First, our kernel represents liquid-topped clouds, while their kernel represents clouds of all phases. Second, our kernel is partitioned by re and LWP, while their kernel is partitioned by cloud-top pressure and cloud optical thickness. Besides these exceptions, we calculate the kernel by closely following the method of Zelinka et al. (2012).

The first step of the kernel calculation is to quantify the clear-sky upward SW flux at the top of the atmosphere for various combinations of surface albedo, latitude, and calendar month. Calculations are performed with the RRTMG radiative transfer model (Clough et al., 2005) using inputs that include the climatological seasonal cycle of humidity from MERRA-2, a standard ozone profile, and a solar constant of 1361 W m−2. For a given latitude and month, we chose a day in the middle of the month and calculate the average of the cosine of the solar zenith angle μi for each 1 h interval throughout the day. We then scale μi by SW,CERES/SW,day, where SW,day is the daily-mean insolation for the day in the middle of the month, and SW,CERES is the monthly-mean insolation from CERES satellite data (Loeb et al., 2018). This step ensures that the monthly-mean insolation for the radiative kernel is equal to that of CERES. We compute the clear-sky SW flux for each of the 24 μi terms and then average the results. The calculations are performed using surface albedo of 0, 0.5, and 1. The final result is a matrix of clear-sky upward SW flux at the top of the atmosphere as a function of surface albedo, latitude, and calendar month.

The next step is similar to the clear-sky calculations except that an overcast and horizontally uniform liquid cloud is introduced in the radiation code. The re and LWP of the cloud are varied, and cloud-top pressure is set to 850 hPa to match the modal value retrieved by MODIS. For re, we use the standard MODIS retrieval algorithm, re,std, and the 3.7 µm retrieval algorithm, re,3.7. Monthly gridded values of re,std and re,3.7 have a correlation coefficient of 0.92 over ocean, but they are generally different because re,3.7 represents conditions closer to the cloud top (Platnick, 2000). We therefore prescribe re inside the cloud as


where τc is the visible optical depth below cloud top, and r̃e,stdmre,3.7+b, where m=1.14 and b=-0.35µm. The coefficients m and b are determined by regressing re,std against re,3.7 using all monthly 1× 1 ocean grid boxes, and the τc=3 threshold is chosen from weighting functions estimated by Platnick (2000). Using these relationships, we calculate the top-of-atmosphere SW flux with different combinations of re,3.7 and LWP that correspond to the bins of the MODIS re–LWP joint histogram. Separate calculations are performed for synthetic clouds at the four edges of each bin, and the results are averaged to get one value of upward SW flux at the top of the atmosphere for each bin. We then subtract the resulting value from the clear-sky upward SW flux to determine the SW cloud radiative effect (CRE). These calculations produce a matrix of SW CRE above an overcast liquid cloud as a function of latitude, surface albedo, calendar month, re, and LWP.

The final step of the calculation is to convert the overcast-sky SW CRE to a cloud radiative kernel, K. Let Crl be the liquid-cloud fraction in effective radius bin r and LWP bin l. The K matrix represents how anomalies of Crl affect R with all non-cloud factors fixed:


Krl is computed by dividing the overcast-sky SW CRE by 100 %. We apply linear interpolation to transform Krl from latitude–surface-albedo space to latitude–longitude space using the climatological seasonal cycle of clear-sky surface albedo from CERES. The final radiative kernel has units of Wm-2%-1 (watts per square meter per percentage of cloud fraction), and it is a function of latitude, longitude, calendar month, re, and LWP (Fig. 1b).

We validate the radiative kernel by comparing R observations from the CERES FluxByCldTyp dataset (Sun et al., 2022) with R estimates computed from the kernel and MODIS re–LWP joint histograms. The kernel estimates are regressed against the CERES observations using data from all monthly 1× 1 ocean grid boxes between 55 S and 55 N from 2003 through 2020. The regression slope is 1.046 ± 0.005 (95 % CI) (Fig. A1). Thus, biases of our radiative kernel and differences between the MODIS and CERES cloud-phase retrieval algorithms cause the kernel-based values of R to overestimate the magnitude of their CERES counterpart by +4.6 ± 0.5 %.

Figure A1Validation of the SW cloud radiative kernel. The vertical axis shows monthly SW flux anomalies induced by liquid-topped clouds estimated with the radiative kernel and MODIS re–LWP joint histograms (Rkernel). The horizontal axis shows monthly SW flux anomalies induced by liquid-topped clouds observed by CERES (RCERES). Data are plotted as a joint histogram compiled from all monthly 1× 1 grid boxes over ocean between 55 S and 55 N from 2003 through 2020. The color scale is logarithmic, and bins with fewer than 100 counts are shaded white for clarity. The black line is the ordinary least-squares regression fit. The regression slope and its 95 % CI are printed in the bottom-right corner.


Appendix B: Decomposing cloud radiative effects

For a given latitude, longitude, and time, the total liquid-cloud-induced SW flux anomaly at the top of the atmosphere, R, is

(B1) R = r = 1 6 l = 1 7 K r l C r l .

We decompose the term on the right side of Eq. (B1) to estimate how much cloud-amount anomalies, re anomalies, and LWP anomalies contribute to R. The decomposition closely follows the method described in Appendix B of Zelinka et al. (2013) except that our radiative kernel and histogram have different dimensions.

First, let Ctot be the total liquid-cloud fraction summed over all histogram bins. We express Crl as

(B2) C r l = C r l C tot C tot + C r l ,

where overbars denote values from the local climatological seasonal cycle. The first term on the right side of Eq. (B2) represents the anomalies of Crl that would occur if Ctot were distributed among the re–LWP bins such that the normalized distribution in the histogram remains the same as the climatology. In other words, this term accounts for a change in total liquid-cloud fraction, holding the proportion of clouds in each histogram bin fixed. The second term on the right side of Eq. (B2) accounts for anomalies of Crl that remain after removing (Crl/Ctot)Ctot. This term represents shifts in the distribution of re and LWP with the total liquid-cloud fraction fixed. By construction, Crl vanishes when it is summed over all histogram bins.

Next, we decompose the radiative kernel into two terms:

(B3) K r l = K 0 + K ^ r l .

Here, K0 is the average of Krl weighted by the climatological cloud fraction,

(B4) K 0 r = 1 6 l = 1 7 C r l C tot K r l ,

and K^rlKrl-K0. With the relationships in Eqs. (B1)–(B4), R can be expressed as


We next decompose K^rl into three components:






R can then be expressed as

(B5) R = K 0 C tot + r = 1 6 K ^ r l = 1 7 C r l + l = 1 7 K ^ l r = 1 6 C r l + r = 1 6 l = 1 7 K ^ res C r l .

The first term on the right side of Eq. (B5) is the SW flux anomaly that would occur if the anomaly of total liquid-cloud fraction were distributed among the re–LWP bins such that the proportion of cloud fraction in each bin is the same as the climatology. This term represents the contribution of cloud-amount anomalies to R. The second term on the right side results from multiplying an effective kernel that accounts for systematic variations in re by the change in cloud fraction at each re bin. This term represents the contribution of re anomalies to R with LWP and total liquid-cloud fraction held fixed. The third term on the right side is similar to the second term except that it represents the contribution of LWP anomalies to R with re and total liquid-cloud fraction held fixed. The final term on the right is the residual of the decomposition. The cloud amount, re, LWP, and residual components of R are denoted by RCF, Rre, RLWP, and Rres, respectively. In the Supplement, we validate the R decomposition using synthetic-data test cases in which Rre and RLWP can be estimated theoretically with the two-stream radiative transfer approximation (Fig. S3). We also verify that Rre and RLWP are similar when different common assumptions are made about cloud vertical structure, including a vertically uniform cloud model, an adiabatic cloud model, and the two-layer cloud model from the kernel calculation.

The final step of the decomposition is to adjust RCF to account for obscuration effects from non-liquid clouds. Because MODIS is a passive instrument, changes in non-liquid clouds can artificially change the retrieved liquid-cloud fraction if they obscure liquid clouds from the satellite view. We control for these obscuration effects by replacing Ctot in Eq. (B5) by [Ctot/(100%-Itot)](100%-Itot), where Itot is the retrieved fraction of non-liquid clouds. This change of variables is adapted from the procedure recommended by Scott et al. (2020).

We also wish to compare the SW and longwave (LW) components of cloud radiative effects. However, our radiative kernel assumes a constant cloud-top pressure, so it cannot accurately assess the LW component (Appendix A). Instead, we analyze observations of the SW and LW radiative effects of liquid-topped clouds from the CERES FluxByCldTyp dataset. We perform cloud-controlling factor analysis with the CERES data to attain estimates of ERFaci in addition to those presented in the main paper. The CERES-based estimate of ERFaci averaged over ocean between 55 S and 55 N is 1.52 ± 0.51 W m−2 for SW radiation and +0.21 ± 0.13 W m−2 for LW radiation. The LW component offsets 14 ± 7 % of the SW component (95 % CIs). These values determine β in the upper-bound estimate of global net ERFaci (Eq. 2).

Appendix C: Uncertainty

Uncertainty in our ocean-average ERFaci estimates arises from uncertainty in the regression coefficients representing R/lns and uncertainty in the model estimates of Δln s. We first quantify the component that is attributable to regression-coefficient uncertainty. For a particular latitude–longitude grid box i, let ϵi represent the half-width of the 95 % confidence interval of the grid-box mean ERFaci. We estimate ϵi as


where σi is the standard error of the regression coefficient, Nnom,i is the nominal number of temporal degrees of freedom (the number of months in the record), Neff,i is the effective number of temporal degrees of freedom, square brackets indicate the central estimate of a parameter, and ti is the critical value of Student's t distribution at the (1-α/2)100% significance level using Neff,i-8 degrees of freedom and α=0.05. The ratio Nnom,i/Neff,i is estimated as (1+a)/(1-a), where a is the temporal lag-1 autocorrelation of Ri. The ϵi terms are then combined to account for spatial averaging over the domain. Uncertainty of the domain-average forcing, δobs, is


where Nnom is the nominal number of spatial degrees of freedom (the number of latitude–longitude grid boxes in the domain), Neff is the effective number of spatial degrees of freedom, and wi is the ocean area in grid box i. The ratio Nnom/Neff is estimated by applying Eq. (5) of Bretherton et al. (1999) to the gridded R data. The resulting value of δobs represents the half-width of the 95 % confidence interval of ERFaci that is attributable to regression-coefficient uncertainty. Confidence intervals for the spatial averages of R/lns and R/lnÑd are calculated similarly.

The second source of uncertainty of ERFaci arises from inter-model spread in the estimates of Δln s. Because we have estimates from 20 climate models, we construct a 95 % confidence interval for ERFaci that excludes 1 model and encompasses the range of the other 19. We first calculate 20 estimates of ERFaci by multiplying Δln s from each model by [R/lns]. The half-width of the confidence interval, δln s, is estimated as the minimum of |c1-c19|/2 and |c2-c20|/2, where c1, c2, c19, and c20 are the smallest, second smallest, second largest, and largest values of the 20 ERFaci estimates, respectively. This uncertainty analysis accounts for inter-model differences in aerosol processing, but it does not account for uncertainty in anthropogenic sulfur-dioxide emissions because the climate model simulations apply the same emission values. However, sulfur-dioxide emissions depend primarily on the sulfur content of fuel rather than the conditions of combustion, so global emission inventories have a relatively small uncertainty of about ± 11 % (90 % CI from Smith et al., 2011). This is expected to cause an equivalent fractional uncertainty in the global burden of anthropogenic sulfate aerosol (Charlson et al., 1992; Stevens, 2015). In contrast, inter-model differences in aerosol processing lead to an uncertainty of ± 43 % in the change in s since 1850 averaged over ocean (90 % CI estimated by computing the interval that includes 18 of the 20 CMIP6 models). The quadrature sum of these two components determines their combined uncertainty, so uncertainty in aerosol processing dominates the overall uncertainty in global-mean Δln s. Approximating δln s from inter-model differences is therefore justified. Finally, we note that ERFaci strongly depends on the preindustrial aerosol state, so our estimates are approximately valid as long as the inter-model spread in preindustrial aerosol concentrations encompasses the true values (Carslaw et al., 2013; Kinne, 2019). The possibility that CMIP6 models have systematic biases that violate this condition cannot be ruled out at this time.

The estimated δobs and δln s represent independent sources of uncertainty, so they are combined in quadrature. The overall 95 % confidence interval is given by ERFaci±δobs2+δlns2. Confidence intervals for IRFaci, ALWP, ACF, and the total cloud adjustment are calculated similarly.

Appendix D: Filtering of partly cloudy pixels

Different filtering methods for partly cloudy pixels in the MODISCLD and MODISCLD+PCL cases lead to offsetting differences in the estimated LWP and cloud-fraction adjustments (Fig. 3). One possible explanation for this discrepancy is that CCN anomalies cause changes in the morphology or horizontal size of liquid clouds, thereby changing the relative amounts of partly and fully cloudy pixels (Possner et al., 2018; Rosenfeld et al., 2006). Here we examine a case study to demonstrate this concept.

We analyze instantaneous pixel data from the MODIS MOD06_L2 dataset collection 6.1 (Platnick et al., 2015) obtained during a single overpass of the Terra satellite on 27 September 2019. On this day, MODIS measured stratocumulus clouds in the southeast Pacific Ocean with different forms of mesoscale cellular convection. The clouds in Box C in Fig. D1a mostly exhibit closed cells, and the clouds in Box O mostly exhibit open cells. We select data from these boxes and bin the liquid-cloud pixels into histograms of cloud fraction partitioned by LWP. Let Cl represent cloud fraction in LWP bin l. In one case, Cl and Ctot are computed by counting only the fully cloudy pixels (CLD), and in a second case, Cl and Ctot are computed by counting both fully and partly cloudy pixels (CLD+PCL). Box C contains mostly fully cloudy pixels, so the two cases have similar values of Cl and Ctot (Fig. D1b). In contrast, Box O contains broken clouds that have a smaller horizontal scale, a larger perimeter-to-area ratio, and a larger proportion of partly cloudy pixels. The partly cloudy pixels cover about 13 % of the area of Box O, and their retrieved LWP is usually smaller than that of the fully cloudy pixels in the box. This causes a difference in Ctot and the Cl distribution between the CLD and CLD + PCL cases (Fig. D1c). Filtering of partly cloudy pixels therefore affects the grid-box-level statistics of cloud fraction and LWP in this example.

We next examine differences between the cloud-fraction histograms of the two boxes to demonstrate the implications for estimating R. For the purpose of this demonstration, let the baseline cloud population be defined by the clouds in Box O, and let the cloud-fraction anomalies be defined by the cloud fraction in Box C minus the cloud fraction in Box O. The baseline and anomalies are indicated by overbars and primes, respectively. Cl can be decomposed according to


where ClCl-ClCtotCtot. This decomposition is equivalent to Eq. (B2) except that it is performed on a one-dimensional LWP histogram rather than a two-dimensional LWP–re joint histogram. The first term on the right side of the equation determines RCF, and the second term determines RLWP. Compared with the CLD case, the CLD + PCL case has a smaller value of Ctot (Fig. D1d), which reduces the magnitude of RCF. Furthermore, because the partly cloudy pixels in Box O occupy the smallest LWP bins, including these pixels in the histogram causes a less extreme shift in the Cl distribution towards small LWP values between Box O and Box C. This reduces the magnitude of RLWP in the CLD + PCL case relative to that of the CLD case. Thus, including partly cloudy pixels in the histograms leads to offsetting changes in RCF and RLWP that reduce the magnitude of both terms. Our estimates of global cloud adjustments depend on filtering of partly cloudy pixels in a similar way, suggesting that the adjustments may involve changes in cloud size or morphology as well.

Figure D1Case study demonstrating how different filtering methods for partly cloud-covered pixels lead to different estimates of LWP and cloud-fraction anomalies when the cloud morphology changes. (a) Visible image of stratocumulus clouds over the southeast Pacific Ocean taken on 27 September 2019 from MODIS on the Terra satellite. Most clouds in Box C exhibit closed mesoscale cellular convection, and most clouds in Box O exhibit open mesoscale cellular convection. Both boxes span 6 latitude and 6 longitude. (b) Liquid-cloud fraction partitioned by LWP (Cl) in Box C. In case CLD, the Cl histogram includes only fully cloud-covered pixels, and in case CLD + PCL, the histogram includes both fully and partly cloud-covered pixels. The CLD and CLD + PCL cases are shown with blue and dashed black lines, respectively. The total liquid-cloud fraction (Ctot) is printed on the figure. (c, d) Similar to (b) except that (c) shows Box O and (d) shows the difference between Box C and Box O. Primes in (d) represent the Box-C average minus the Box-O average.

Code availability

MATLAB code used to analyze data can be obtained by contacting the corresponding author. The RRTMG radiative transfer model is publicly available at (Atmospheric and Environmental Research, 2020).

Data availability

All satellite data, reanalysis, and GCM output used in this study are publicly available. MODIS cloud histograms are available from the National Aeronautics and Space Administration (NASA) Level-1 and Atmosphere Archive and Distribution System (, NASA, 2023b), MODIS Nd data are available from the Centre for Environmental Data Analysis (, Gryspeerdt et al., 2022b), CERES data are available from the CERES online ordering tool ​(, CERES, 2023a;, CERES, 2023b), MERRA-2 data are available from the NASA Goddard Earth Sciences Data and Information Services Center (, NASA, 2023a), and CMIP6 output is available from the CMIP6 online archive (, WCRP, 2022). The SW cloud radiative kernel is available from GitHub (2023,, Wall, 2023a;, Wall, 2023b). Estimates of the components of aerosol indirect effects and ERFaci averaged over ocean are listed in Table S3.


The supplement related to this article is available online at:

Author contributions

CJW designed research, performed research, analyzed data, and wrote the original paper draft. TS and AP contributed ideas that shaped the study. All authors helped revise the original paper draft.

Competing interests

The contact author has declared that none of the authors has any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank Daniel McCoy, Ed Gryspeerdt, Robert Pincus, and Velle Toll for helpful discussions; Paul Hubanks, Robert Pincus, and Steven Platnick for supporting the design and implementation of the MODIS joint histograms; and Ying Chen and Jianhao Zhang for providing constructive reviews that improved the manuscript.

Financial support

This project has received funding from the European Union’s Horizon 2020 Research and Innovation program under the Marie Skłodowska-Curie Action (grant agreement no. 101019911).

Review statement

This paper was edited by Matthew Lebsock and reviewed by Ying Chen and Jianhao Zhang.


Albrecht, B. A.: Aerosols, cloud microphysics, and fractional cloudiness, Science, 245, 1227–1230,, 1989. 

Andersen, H., Cermak, J., Fuchs, J., Knutti, R., and Lohmann, U.: Understanding the drivers of marine liquid-water cloud occurrence and properties with global observations using neural networks, Atmos. Chem. Phys., 17, 9535–9546,, 2017. 

Andreae, M. O. and Rosenfeld, D.: Aerosol-Cloud-Precipitation Interactions Part 1: The Nature and Sources of Cloud-Active Aerosols, Earth-Sci. Rev., 89, 13–41,, 2008. 

Atmospheric and Environmental Research: RRTMG-SW, GitHub [code], (last access: 11 March 2021), 2020. 

Bellouin, N., Quaas, J., Gryspeerdt, E., Kinne, S., Stier, P., Watson-Parris, D., Boucher, O., Carslaw, K. S., Christensen, M., Daniau, A.-L., Dufresne, J.-L., Feingold, G., Fiedler, S., Forster, P., Gettelman, A., Haywood, J. M., Lohmann, U., Malavelle, F., Mauritsen, T., McCoy, D. T., Myhre, G., Mülmenstädt, J., Neubauer, D., Possner, A., Rugenstein, M., Sato, Y., Schulz, M., Schwartz, S. E., Sourdeval, O., Storelvmo, T., Toll, V., Winker, D., and Stevens, B.: Bounding Global Aerosol Radiative Forcing of Climate Change, Rev. Geophys., 58, e2019RG000660,, 2020. 

Bender, F. A. M., Frey, L., McCoy, D. T., Grosvenor, D. P., and Mohrmann, J. K.: Assessment of aerosol–cloud–radiation correlations in satellite observations, climate models and reanalysis, Clim. Dynam., 52, 4371–4392,, 2019. 

Bretherton, C. S., Widmann, M., Dymnikov, V. P., Wallace, J. M., and Bladé, I.: Effective Number of Degrees of Freedom of a Spatial Field, J. Climate, 12, 1990–2009<1990:TENOSD>2.0.CO;2, 1999. 

Bretherton, C. S., Blossey, P. N., and Uchida, J.: Cloud Droplet Sedimentation, Entrainment Efficiency, and Subtropical Stratocumulus Albedo, Geophys. Res. Lett., 34, 1–5,, 2007. 

Carslaw, K. S., Lee, L. A., Reddington, C. L., Pringle, K. J., Rap, A., Forster, P. M., Mann, G. W., Spracklen, D. V., Woodhouse, M. T., Regayre, L. A., and Pierce, J. R.: Large contribution of natural aerosols to uncertainty in indirect forcing, Nature, 503, 67–71,, 2013. 

Charlson, R. J., Schwartz, S. E., Hales, J. M., Cess, R. D., Coakley, J. A., Hansen, J. E., and Hofmann, D. J.: Climate Forcing by Anthropogenic Aerosols, Science, 255, 423–430,, 1992. 

CERES: EBAF_Ed4.1, NASA [data set], (last access: 1 January 2023), 2023a. 

CERES: FluxByCldTyp_Ed4.1, NASA [data set], (last access: 1 January 2023), 2023b. 

Chen, Y., Haywood, J., Wang, Y., Malavelle, F., Jordan, G., Partridge, D., Fieldsend, J., De Leeuw, J., Schmidt, A., Cho, N., Oreopoulos, L., Platnick, S., Grosvenor, D., Field, P., and Lohmann, U.: Machine Learning Reveals Climate Forcing from Aerosols is Dominated by Increased Cloud Cover, Nat. Geosci., 15, 609–614,, 2022. 

Chen, Y. C., Christensen, M. W., Stephens, G. L., and Seinfeld, J. H.: Satellite-Based Estimate of Global Aerosol-Cloud Radiative Forcing by Marine Warm Clouds, Nat. Geosci., 7, 643–646,, 2014. 

Clough, S. A., Shephard, M. W., Mlawer, E. J., Delamere, J. S., Iacono, M. J., Cady-Pereira, K., Boukabara, S., and Brown, P. D.: Atmospheric Radiative Transfer Modeling: A Summary of the AER Codes, J. Quant. Spectrosc. Ra., 91, 233–244,, 2005. 

Dagan, G., Koren, I., Altaratz, O., and Heiblum, R. H.: Time-dependent, non-monotonic response of warm convective cloud fields to changes in aerosol loading, Atmos. Chem. Phys., 17, 7435–7444,, 2017. 

Diamond, M. S., Director, H. M., Eastman, R., Possner, A., and Wood, R.: Substantial Cloud Brightening From Shipping in Subtropical Low Clouds, AGU Adv., 1, 1–28,, 2020. 

Feingold, G., Goren, T., and Yamaguchi, T.: Quantifying albedo susceptibility biases in shallow clouds, Atmos. Chem. Phys., 22, 3303–3319,, 2022. 

Forster, P., Storelvmo, T., Armour, K., Collins, W., Dufresne, J.-L., Frame, D., Lunt, D. J., Mauritsen, T., Palmer, M. D., Watanabe, M., Wild, M., and Zhang, H.: The Earth's Energy Budget, Climate Feedbacks, and Climate Sensitivity, in: Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S. L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis, M. I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou, B., Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA,, pp. 923–1054, 2021. 

Gelaro, R., McCarty, W., Suárez, M. J., Todling, R., Molod, A., Takacs, L., Randles, C. A., Darmenov, A., Bosilovich, M. G., Reichle, R., Wargan, K., Coy, L. Cullather, R., Draper, C., Akella, S., Buchard, V., Conaty, A., da Silva, A. M., Gu, W., Kim, G. K., Koster, R., Lucchesi, R., Merkova, D., Nielsen, J. E., Partyka, G., Pawson, S., Putman, W., Rienecker, M., Schubert, S. D., Sienkiewicz, M., and Zhao, B.: The Modern-Era Retrospective Analysis for Research and Applications, version 2 (MERRA-2), J. Climate, 30, 5419–5454,, 2017. 

Ghan, S., Wang, M., Zhang, S., Ferrachat, S., Gettelman, A., Griesfeller, J., Kipling, Z., Lohmann, U., Morrison, H., Neubauer, D., Partridge, D. G., Stier, P., Takemura, T., Wang, H., and Zhang, K.: Challenges in Constraining Anthropogenic Aerosol Effects on Cloud Radiative Forcing using Present-Day Spatiotemporal Variability, P. Natl. Acad. Sci. USA, 113, 5804–5811,, 2016. 

Glassmeier, F., Hoffmann, F., Johnson, J. S., Yamaguchi, T., Carslaw, K. S., and Feingold, G.: Aerosol-Cloud-Climate Cooling Overestimated by Ship-Track Data, Science, 371, 485–489,, 2021. 

Goren, T., Feingold, G., Gryspeerdt, E., Kazil, J., Kretzschmar, J., Jia, H., and Quaas, J.: Projecting Stratocumulus Transitions on the Albedo – Cloud Fraction Relationship Reveals Linearity of Albedo to Droplet Concentrations, Geophys. Res. Lett., 49, e2022GL101169,, 2022. 

Goren, T., Sourdeval, O., Kretzschmar, J., and Quaas, J.: Spatial Aggregation of Satellite Observations Leads to an Overestimation of the Radiative Forcing Due To Aerosol-Cloud Interactions, Geophy. Res. Lett., 50, e2023GL105282,, 2023. 

Grosvenor, D. P., Sourdeval, O., Zuidema, P., Ackerman, A., Alexandrov, M. D., Bennartz, R., Boers, R., Cairns, B., Chiu, J. C., Christensen, M., Deneke, H., Diamond, M., Feingold, G., Fridlind, A., Hünerbein, A., Knist, C., Kollias, P., Marshak, A., McCoy, D. T, Merk, D., Painemal, D., Rausch, J., Rosenfeld, D., Russchenberg, H., Seifert, P., Sinclair, K., Stier, P., van Diedenhoven, B., Wendisch, M., Werner, F., Wood, R., Zhang, Z., and Quaas, J.: Remote Sensing of Droplet Number Concentration in Warm Clouds: A Review of the Current State of Knowledge and Perspectives, Rev. Geophys., 56, 409–453,, 2018. 

Gryspeerdt, E., Quaas, J., and Bellouin, N.: Constraining the Aerosol Influence on Cloud Fraction, J. Geophys. Res.-Atmos., 121, 3566–3583,, 2016. 

Gryspeerdt, E., Goren, T., Sourdeval, O., Quaas, J., Mülmenstädt, J., Dipu, S., Unglaub, C., Gettelman, A., and Christensen, M.: Constraining the aerosol influence on cloud liquid water path, Atmos. Chem. Phys., 19, 5331–5347,, 2019. 

Gryspeerdt, E., Mülmenstädt, J., Gettelman, A., Malavelle, F. F., Morrison, H., Neubauer, D., Partridge, D. G., Stier, P., Takemura, T., Wang, H., Wang, M., and Zhang, K.: Surprising similarities in model and observational aerosol radiative forcing estimates, Atmos. Chem. Phys., 20, 613–623,, 2020. 

Gryspeerdt, E., McCoy, D. T., Crosbie, E., Moore, R. H., Nott, G. J., Painemal, D., Small-Griswold, J., Sorooshian, A., and Ziemba, L.: The impact of sampling strategy on the cloud droplet number concentration estimated from satellite data, Atmos. Meas. Tech., 15, 3875–3892,, 2022a. 

Gryspeerdt, E., McCoy, D., Crosbie, E., Moore, R. H., Nott, G. J., Painemal, D., Small-Griswold, J., Sorooshian, A., and Ziemba, L.: Cloud droplet number concentration, calculated from the MODIS cloud optical properties retrieval and gridded using different sampling strategies, NERC EDS Centre for Environmental Data Analysis [data set],, 2022b. 

IPCC: Annex III: Tables of historical and projected well-mixed greenhouse gas mixing ratios and effective radiative forcing of all climate forcers, edited by: Dentener F. J., Hall, B., and Smith, C., in: Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, Masson-Delmotte, V., Zhai, P., Pirani, A., Connors, S. L., Péan, C., Berger, S., Caud, N., Chen, Y., Goldfarb, L., Gomis, M. I., Huang, M., Leitzell, K., Lonnoy, E., Matthews, J. B. R., Maycock, T. K., Waterfield, T., Yelekçi, O., Yu, R., and Zhou , B., Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA,, pp. 2139–2152, 2021. 

Jia, H., Ma, X., Yu, F., and Quaas, J.: Significant Underestimation of Radiative Forcing by Aerosol–Cloud Interactions Derived from Satellite-Based Methods, Nat. Commun., 12, 1–11,, 2021. 

Kang, L., Marchand, R. T., Wood, R., and McCoy, I. L.: Coalescence Scavenging Drives Droplet Number Concentration in Southern Ocean Low Clouds, Geophys. Res. Lett., 49, 1–10,, 2022. 

Kinne, S.: Aerosol radiative effects with MACv2, Atmos. Chem. Phys., 19, 10919–10959,, 2019. 

Kramer, R. J., He, H., Soden, B. J., Oreopoulos, L., Myhre, G., Forster, P. M., and Smith, C. J.: Observational Evidence of Increasing Global Radiative Forcing, Geophys. Res. Lett., 48, 1–11,, 2021. 

Lewis, H., Bellon, G., and Dinh, T.: Upstream Large-Scale Control of Subtropical Low-Cloud Climatology, J. Climate, 36, 3289–3303,, 2023. 

Loeb, N. G., Doelling, D. R., Wang, H., Su, W., Nguyen, C., Corbett, J. G., Liang, L., Mitrescu, C., Rose, F. G., and Kato, S.: Clouds and the Earth's Radiant Energy System (CERES) Energy Balanced and Filled (EBAF) Top-of-Atmosphere (TOA) Edition-4.0 Data Product, J. Climate, 31, 895–918,, 2018. 

Ma, P., Rasch, P. J., Chepfer, H., Winker, D. M., and Ghan, S. J.: Observational Constraint on Cloud Susceptibility Weakened by Aerosol Retrieval Limitations, Nat. Commun., 9, 1–10,, 2018. 

Manshausen, P., Watson-Parris, D., Christensen, M. W., Jalkanen, J. P., and Stier, P.: Invisible ship tracks show large cloud sensitivity to aerosol, Nature, 610, 101–106,, 2022. 

Matthews, H. D. and Zickfeld, K.: Climate Response to Zeroed Emissions of Greenhouse Gases and Aerosols, Nat. Clim. Change, 2, 338–341,, 2012. 

Mauger, G. S. and Norris, J. R.: Meteorological Bias in Satellite Estimates of Aerosol-Cloud Relationships, Geophys. Res. Lett., 34, 1–5,, 2007. 

Mauritsen, T. and Pincus, R.: Committed Warming Inferred from Observations, Nat. Clim. Change, 7, 652–655,, 2017. 

McCoy, D. T., Bender, F. A. M., Mohrmann, J. K. C., Hartmann, D. L., Wood, R., and Grosvenor, D. P.: The global aerosol-cloud first indirect effect estimated using MODIS, MERRA, and AeroCom. J. Geophys. Res.- Atmos., 122, 1779–1796,, 2017. 

McCoy, D. T., Field, P. R., Schmidt, A., Grosvenor, D. P., Bender, F. A.-M., Shipway, B. J., Hill, A. A., Wilkinson, J. M., and Elsaesser, G. S.: Aerosol midlatitude cyclone indirect effects in observations and high-resolution simulations, Atmos. Chem. Phys., 18, 5821–5846,, 2018. 

Myers, T. A., Scott, R. C., Zelinka, M. D., Klein, S. A., Norris, J. R., and Caldwell, P. M.: Observational Constraints on Low Cloud Feedback Reduce Uncertainty of Climate Sensitivity, Nat. Clim. Change, 11, 501–507,, 2021. 

NASA: MERRA-2, Goddard Earth Sciences Data and Information Services Center [data set],, 2023a. 

NASA: MCD06COSP_M3_MODIS – MODIS (Aqua/Terra) Cloud Properties Level 3 monthly, 1x1 degree grid, Distributed Active Archive Center [data set],, 2023b. 

Naud, C. M., Posselt, D. J., and van den Heever, S. C.: Observed Covariations of Aerosol Optical Depth and Cloud Cover in Extratropical Cyclones, J. Geophys. Res.-Atmos., 122, 10338–10356,, 2017. 

Painemal, D., Chiu, J.-Y. C., Minnis, P., Yost, C., Zhou, X., Cadeddu, M., Eloranta, E., Lewis, E. R., Ferrare, R., and Kollias, P.: Aerosol and Cloud Microphysics Covariability in the Northeast Pacific Boundary Layer Estimated with Ship-Based and Satellite Remote Sensing Observations, J. Geophys. Res.-Atmos., 122, 2403–2418,, 2017. 

Painemal, D. and Zuidema, P.: Assessment of MODIS Cloud Effective Radius and Optical Thickness Retrievals over the Southeast Pacific with VOCALS-REx In Situ Measurements, J. Geophys. Res.-Atmos., 116, 1–16,, 2011. 

Pincus, R. and Baker, M. B.: Effect of Precipitation on the Albedo Susceptibility of Clouds in the Marine Boundary Layer, Nature, 372, 250–252,, 1994. 

Pincus, R., Hubanks, P. A., Platnick, S., Meyer, K., Holz, R. E., Botambekov, D., and Wall, C. J.: Updated observations of clouds by MODIS for global model assessment, Earth Syst. Sci. Data, 15, 2483–2497,, 2023. 

Platnick, S.: Vertical Photon Transport in Cloud Remote Sensing Problems, J. Geophys. Res.-Atmos., 105, 22919–22935,, 2000. 

Platnick, S., Ackerman, S., and King, M.: MODIS Atmosphere L2 Cloud Product (06_L2). NASA MODIS Adaptive Processing System, Goddard Space Flight Center, USA,, 2015. 

Platnick, S., Meyer, K. G., King, M. D., Wind, G., Amarasinghe, N., Marchant, B., Arnold, G. T., Zhang, Z., Hubanks, P. A., Holz, R. E., Yang, P., Ridgway, W. L., and Reidi, J.: The MODIS Cloud Optical and Microphysical Products: Collection 6 Updates and Examples from Terra and Aqua, IEEE T. Geosci. Remote, 3, 163–186,, 2017. 

Possner, A., Wang, H., Wood, R., Caldeira, K., and Ackerman, T. P.: The efficacy of aerosol–cloud radiative perturbations from near-surface emissions in deep open-cell stratocumuli, Atmos. Chem. Phys., 18, 17475–17488,, 2018. 

Possner, A., Eastman, R., Bender, F., and Glassmeier, F.: Deconvolution of boundary layer depth and aerosol constraints on cloud water path in subtropical stratocumulus decks, Atmos. Chem. Phys., 20, 3609–3621,, 2020. 

Raghuraman, S. P., Paynter, D., and Ramaswamy, V.: Anthropogenic forcing and response yield observed positive trend in Earth's energy imbalance, Nat. Commun., 12, 1–10,, 2021. 

Randles, C. A., da Silva, A. M., Buchard, V., Colarco, P. R., Darmenov, A., Govindaraju, R., Smirnov, A., Holben, B., Ferrare, R., Hair, J., Shinozuka, Y., and Flynn, C. J.: The MERRA-2 Aerosol Reanalysis, 1980 Onward. Part I: System Description and Data Assimilation Evaluation, J. Climate, 30, 6823–6850,, 2017. 

Rosenfeld, D., Kaufman, Y. J., and Koren, I.: Switching cloud cover and dynamical regimes from open to closed Benard cells in response to the suppression of precipitation by aerosols, Atmos. Chem. Phys., 6, 2503–2511,, 2006. 

Rosenfeld, D., Andreae, M. O., Asmi, A., Chin, M., Leeuw, G., Donovan, D. P., Kahn, R., Kinne, S., Kivekäs, N., Kulmala, M., Lau, W., Schmidt, S. K., Suni, T., Wagner, T., Wild, M., and Quaas, J.: Global Observations of Aerosol-Cloud-Precipitation-Climate Interactions. Rev. Geophys., 52, 750–808,, 2014. 

Rosenfeld, D., Zhu, Y., Wang, M., Zheng, Y., Goren, T., and Yu, S.: Aerosol-Driven Droplet Concentrations Dominate Coverage and Water of Oceanic Low-Level Clouds, Science, 363, eaav0566,, 2019. 

Scott, R. C., Myers, T. A., Norris, J. R., Zelinka, M. D., Klein, S. A., Sun, M., and Doelling, D. R.: Observed Sensitivity of Low-Cloud Radiative Effects to Meteorological Perturbations over the Global Oceans, J. Climate, 33, 7717–7734,, 2020. 

Seifert, A., Heus, T., Pincus, R., and Stevens, B.: Large-Eddy Simulation of the Transient and Near-Equilibrium Behavior of Precipitating Shallow Convection, J. Adv. Model. Earth Sy., 6, 513–526,, 2015. 

Sherwood, S. C., Webb, M. J., Annan, J. D., Armour, K. C., Forster, P. M., Hargreaves, J. C., Hegerl, G., Klein, S. A., Marvel, K. D., Rohling, E. J., Watanabe, M., Andrews, T., Braconnot, P., Bretherton, C. S., Foster, G. L., Hausfather, Z., von der Heydt, A. S., Knutti, R., Mauritsen, T., Norris, J. R., Proistosescu, C., Rugenstein, M., Schmidt, G. A., Tokarska, K. B., and Zelinka, M. D.: An Assessment of Earth's Climate Sensitivity Using Multiple Lines of Evidence. Rev. Geophys., 58, 1–92,, 2020. 

Small, J. D., Chuang, P. Y., Feingold, G., and Jiang, H.: Can Aerosol Decrease Cloud Lifetime? Geophys. Res. Lett., 36, 1–5,, 2009. 

Smith, C. J., Harris, G. R., Palmer, M. D., Bellouin, N., Collins, W., Myhre, G., Schulz, M., Golaz, J. C., Ringer, M., Storelvmo, T., and Forster, P. M.: Energy Budget Constraints on the Time History of Aerosol Forcing and Climate Sensitivity, J. Geophys. Res.-Atmos., 126, 1–20,, 2021. 

Smith, S. J., van Aardenne, J., Klimont, Z., Andres, R. J., Volke, A., and Delgado Arias, S.: Anthropogenic sulfur dioxide emissions: 1850–2005, Atmos. Chem. Phys., 11, 1101–1116,, 2011. 

Stevens, B.: Rethinking the Lower Bound on Aerosol Radiative Forcing, J. Climate, 28, 4794–4819,, 2015. 

Sun, M., Doelling, D. R., Loeb, N. G., Scott, R. C., Wilkins, J., Nguyen, L. T., and Mlynczak, P.: Clouds and the Earth's Radiant Energy System (CERES) FluxByCldTyp Edition 4 Data Product, J. Atmos. Ocean. Tech., 39, 303–318,, 2022. 

Toll, V., Christensen, M., Quaas, J., and Bellouin, N.: Weak average liquid-cloud-water response to anthropogenic aerosols, Nature, 572, 51–55,, 2019. 

Twomey, S.: Influence of Pollution on Shortwave Albedo of Clouds, J. Atmos. Sci., 34, 1149–1152,<1149:TIOPOT>2.0.CO;2, 1977. 

Wall, C. J.: MODIS Re-LWP cloud radiative kernel, GitHub [data set], (last access: 10 October 2023), 2023a. 

Wall, C. J.: caseywall7926/MODIS_Re-LWP_kernel: MODIS Re-LWP cloud radiative kernel (v1.0), Zenodo,, 2023b. 

Wall, C. J., Norris, J. R., Possner, A., McCoy, D. T., McCoy, I. L., and Lutsko, N. J.: Assessing Effective Radiative Forcing from Aerosol–Cloud Interactions over the Global Ocean, P. Natl. Acad. Sci. USA, 119,, 2022. 

Watson-Parris, D. and Smith, C. J.: Large uncertainty in future warming due to aerosol forcing, Nat. Clim. Change, 12, 1111–1113,, 2022. 

WCRP: CMIP6 historical simulations, CMIP6 online archive [data set], (last access: 1 May 2022), 2022. 

Wilks, D. S.: “The Stippling Shows Statistically Significant Grid Points.” How Research Results are Routinely Overstated and Overinterpreted, and What to Do about It, B. Am. Meteorol. Soc., 97, 2263–2274,, 2016. 

Wood, R. and Bretherton, C. S.: On the Relationship Between Stratiform Low Cloud Cover and Lower-Tropospheric Stability, J. Climate, 19, 6425–6432,, 2006. 

Wood, R., Leon, D., Lebsock, M., Snider, J., and Clarke, A. D.: Precipitation Driving of Droplet Concentration Variability in Marine Low Clouds, J. Geophys. Res.-Atmos., 117, 1–11,, 2012. 

Zelinka, M. D., Klein, S. A., and Hartmann, D. L.: Computing and Partitioning Cloud Feedbacks Using Cloud Property Histograms. Part I: Cloud Radiative Kernels, J. Climate, 25, 3715–3735,, 2012. 

Zelinka, M. D., Klein, S. A., Taylor, K. E., Andrews, T., Webb, M. J., Gregory, J. M., and Forster, P. M.: Contributions of Different Cloud Types to Feedbacks and Rapid Adjustments in CMIP5, J. Climate, 26, 5007–5027,, 2013. 

Zelinka, M. D., Andrews, T., Forster, P. M., and Taylor, K. E.: Quantifying Components of Aerosol-Cloud-Radiation Interactions in Climate Models, J. Geophys. Res.-Atmos., 119, 7599–7615,, 2014. 

Zhang, J. and Feingold, G.: Distinct regional meteorological influences on low-cloud albedo susceptibility over global marine stratocumulus regions, Atmos. Chem. Phys., 23, 1073–1090,, 2023. 

Zhang, J., Zhou, X., Goren, T., and Feingold, G.: Albedo susceptibility of northeastern Pacific stratocumulus: the role of covarying meteorological conditions, Atmos. Chem. Phys., 22, 861–880,, 2022. 

Zhang, Z. and Platnick, S.: An Assessment of Differences Between Cloud Effective Particle Radius Retrievals for Marine Water Clouds from three MODIS Spectral Bands, J. Geophys. Res.-Atmos., 116, 1–19,, 2011. 

Zhou, X., Zhang, J., and Feingold, G.: On the Importance of Sea Surface Temperature for Aerosol-Induced Brightening of Marine Clouds and Implications for Cloud Feedback in a Future Warmer Climate, Geophys. Res. Lett., 48,, 2021. 

Zhu, Y., Rosenfeld, D., and Li, Z.: Under What Conditions Can We Trust Retrieved Cloud Drop Concentrations in Broken Marine Stratocumulus? J. Geophys. Res. Atmos., 123, 8754–8767,, 2018. 

Executive editor
One of the largest sources of uncertainty in the overall anthropogenic forcing of climate is still the aerosol impact on liquid clouds. Disentangling the various aerosol-cloud interactions helps to improve estimates of the magnitude of global warming in the future. The current study provides the most rigorous method to date in assessing the aerosol radiative effects from satellite observations across the global ocean. The aerosol responses are decomposed into the Twomey effect (cooling due to an increase in cloud-droplet number concentration), and the adjustments of the cloud liquid water path and cloud fraction (often analysed separately) at a near-global scale. The total effective radiative forcing of liquid clouds since 1850 has been found to be negative, with the cloud adjustments larger than the Twomey effect, which was previously thought to be larger.
Short summary
Interactions between aerosol pollution and liquid clouds are one of the largest sources of uncertainty in the effective radiative forcing of climate over the industrial era. We use global satellite observations to decompose the forcing into components from changes in cloud-droplet number concentration, cloud water content, and cloud amount. Our results reduce uncertainty in these forcing components and clarify their relative importance.
Final-revised paper