Articles | Volume 20, issue 24
Research article
17 Dec 2020
Research article |  | 17 Dec 2020

The decomposition of cloud–aerosol forcing in the UK Earth System Model (UKESM1)

Daniel P. Grosvenor and Kenneth S. Carslaw

Climate variability in the North Atlantic influences processes such as hurricane activity and droughts. Global model simulations have identified aerosol–cloud interactions (ACIs) as an important driver of sea surface temperature variability via surface aerosol forcing. However, ACIs are a major cause of uncertainty in climate forcing; therefore, caution is needed in interpreting the results from coarse-resolution, highly parameterized global models.

Here, we separate and quantify the components of the surface shortwave effective radiative forcing (ERF) due to aerosol in the atmosphere-only version of the UK Earth System Model (UKESM1) and evaluate the cloud properties and their radiative effects against observations. We focus on a northern region of the North Atlantic (NA) where stratocumulus clouds dominate (denoted the northern NA region) and a southern region where trade cumulus and broken stratocumulus dominate (southern NA region). Aerosol forcing was diagnosed using a pair of simulations in which the meteorology is approximately fixed via nudging to analysis; one simulation has pre-industrial (PI) and one has present-day (PD) aerosol emissions. This model does not include aerosol effects within the convective parameterization (but aerosol does affect the clouds associated with detrainment) and so it should be noted that the representation of aerosol forcing for convection is incomplete.

Contributions to the surface ERF from changes in cloud fraction (fc), in-cloud liquid water path (LWPic) and droplet number concentration (Nd) were quantified. Over the northern NA region, increases in Nd and LWPic dominate the forcing. This is likely because the already-high fc there reduces the chances of further large increases in fc and allows cloud brightening to act over a larger region. Over the southern NA region, increases in fc dominate due to the suppression of rain by the additional aerosols. Aerosol-driven increases in macrophysical cloud properties (LWPic and fc) will rely on the response of the boundary layer parameterization, along with input from the cloud microphysics scheme, which are highly uncertain processes.

Model grid boxes with low-altitude clouds present in both the PI and PD dominate the forcing in both regions. In the northern NA, the brightening of completely overcast low cloud scenes (100 % cloud cover, likely stratocumulus) contributes the most, whereas in the southern NA the creation of clouds with fc of around 20 % from clear skies in the PI was the largest single contributor, suggesting that trade cumulus clouds are created in response to increases in aerosol. The creation of near-overcast clouds was also important there.

The correct spatial pattern, coverage and properties of clouds are important for determining the magnitude of aerosol forcing, so we also assess the realism of the modelled PD clouds against satellite observations. We find that the model reproduces the spatial pattern of all the observed cloud variables well but that there are biases. The shortwave top-of-the-atmosphere (SWTOA) flux is overestimated by 5.8 % in the northern NA region and 1.7 % in the southern NA, which we attribute mainly to positive biases in low-altitude fc. Nd is too low by 20.6 % in the northern NA and too high by 21.5 % in the southern NA but does not contribute greatly to the main SWTOA biases. Cloudy-sky liquid water path mainly shows biases north of Scandinavia that reach between 50 % and 100 % and dominate the SWTOA bias in that region.

The large contribution to aerosol forcing in the UKESM1 model from highly uncertain macrophysical adjustments suggests that further targeted observations are needed to assess rain formation processes, how they depend on aerosols and the model response to precipitation in order to reduce uncertainty in climate projections.

1 Introduction

Uncertainty in the radiative forcing (RF) from aerosols is the largest of the climate RF uncertainties over the industrial period (Boucher et al.2013). General circulation models (GCMs) that simulate a small magnitude of cooling from aerosols are able to reproduce the observed temperature record if they have a low climate sensitivity and vice versa, which results in large uncertainties in climate sensitivity and therefore also in temperature change predictions (Andreae et al.2005; Golaz et al.2013). Mülmenstädt and Feingold (2018) suggests that one reason for the lack of progress in reducing the uncertainty in aerosol forcing despite years of research is that there has been a lack of studies that target (via evaluation and improvement) the individual processes that cause aerosol–cloud interactions within GCMs.

Aerosol effective radiative forcing (ERF, which differs from RF in that all physical variables are allowed to respond to perturbations except for those concerning the ocean and sea ice; e.g. see Myhre et al.2013) can be separated into a component due to aerosol radiative interactions (ARIs) that occur outside of clouds (sometimes also known as the direct effect) and a component due to aerosol–cloud interactions (ACIs, or indirect effects). The ACI ERF is often also broken down into two further components. The first is due to an increase in cloud droplet concentration (Nd) alone at constant liquid water content (LWC) and constant cloud fraction (fc), which causes a decrease in the cloud droplet effective radius (re). Here, this is termed ERFNd (or often the Twomey effect; Twomey1977). The second ERF component concerns rapid adjustments of LWC (or the vertical integral of this, which is the liquid water path, LWP) for only the cloudy parts of model grid boxes (termed in-cloud LWP, or LWPic here) and/or adjustments in fc that occur in response to the initial decrease in droplet size associated with the Nd increase. These are termed ERFLWPic and ERFfc, respectively.

The mechanisms that cause the adjustments involve several microphysical and thermodynamical processes (Albrecht1989; Stevens et al.1998; Ackerman et al.2004; Bretherton et al.2007; Hill et al.2009; Berner et al.2013; Feingold et al.2015). Simulation of adjustments in GCMs therefore requires the parameterization of many subgrid-scale processes, which are likely to be more difficult for GCMs to get right than ERFNd, where only the change in re needs to be parameterized. It is therefore desirable to separate ERFNd and the adjustment effects within GCMs so that they can be evaluated (against observations and high-resolution models) and improved individually. The observational constraint of ERFNd (which is likely easier than the constraint of adjustments) might then reduce the overall forcing uncertainty and highlight issues with the adjustment part of the forcing. Separating the adjustments into ERFLWPic and ERFfc components will further facilitate more detailed process-level improvements. A further reason to separate the different components is that current models simulate the same forcing with different combinations of ERFNd and adjustment components (Gryspeerdt et al.2020). Therefore, one aim of this study is to separate and quantify the ERF contributions from Nd, LWPic and fc changes in a GCM. A second aim is to also quantify the amount of aerosol forcing from a GCM that originates from different cloud types and changes between cloud types. This too will allow a more targeted model evaluation and improvement in future work.

The aerosol forcing of GCMs is also important regionally, for example, in the North Atlantic (NA) region, which is the focus of this paper. It has been suggested (Booth et al.2012, hereafter B12) that surface radiative aerosol forcing is the dominant driver of the variability in multi-decadal NA sea surface temperatures (SSTs) for the ocean–atmosphere coupled GCM (the UK Met Office HadGEM2-ES model) that was used in the fifth Coupled Model Intercomparison Project (CMIP5). NA SST variability has been linked to impacts on important climate phenomena such as hurricane and tropical storm activity (Zhang and Delworth2006; Smith et al.2010; Dunstone et al.2013); rainfall anomalies in Europe and North America (Sutton and Hodson2005; Sutton and Dong2012); droughts in the African Sahel and Amazonian regions (Hoerling et al.2006; Knight et al.2006; Ackerley et al.2011); Greenland ice-sheet melt rates (Holland et al.2008; Hanna et al.2012); anomalies in sea levels (McCarthy et al.2015); and the midlatitude jet strength (Woollings et al.2015). For a review of changes in the North Atlantic climate system (with a focus on more recent changes), see Robson et al. (2018).

B12 showed that HadGEM2-ES, which represented aerosol–cloud interactions, reproduced the observed NA SSTs with good fidelity in contrast to the CMIP3 models that mostly did not include aerosol–cloud effects. Furthermore, tests using constant aerosols clearly showed the impact of aerosols upon NA SST variability in HadGEM2-ES. Aerosol ARI forcing was found to be negligible in the NA compared to aerosol indirect forcing. The spatial patterns of the indirect forcing, the downwelling surface shortwave (SW) radiation anomalies and the SST anomalies were all consistent, indicating a link between the three. Moreover, a simulation with fixed SSTs showed similar incoming surface SW to that in the coupled model. This suggested that the simulated surface SW anomalies were not brought about by the modification of SSTs as a result of ocean dynamics or other processes. Thus, the implication from the HadGEM2-ES model is that aerosol indirect forcing has a direct local impact on SSTs via the modification of surface downwelling SW radiation.

The claims made in B12 are based upon a GCM and not direct observations. As mentioned earlier, the aerosol forcing in GCMs is highly uncertain for many potential reasons. For example, B12 used a coarse N96 model resolution, and thus the model relies upon parameterizations to represent subgrid cloud formation and cloud–aerosol interactions. It could be the case that HadGEM2-ES overestimates the magnitude of the aerosol forcing and thus overstates its role in driving NA SST variability. Zhang et al. (2013) suggest that the HadGEM2-ES model has some shortcomings in its representation of the ocean that may affect its ability to properly simulate the influence of the ocean. Other papers also argue for an important role for ocean processes in determining the NA SST variability (Ba et al.2014; Knight2005; Menary et al.2015; Robson et al.2016; R. Zhang et al.2016; Yan et al.2018), which may indicate a lesser role for aerosol than simulated in B12.

The aim of this paper is to quantify the mechanisms by which a global climate model (an improved version of the model used in B12) produces aerosol forcing, with a focus on the NA region. Some work breaking down the aerosol ERF of a different GCM into that due to changes in Nd, LWPic and fc has already been recently published (Mülmenstädt et al.2019). They found that in the ECHAM-HAMMOZ model ERFLWPic was quite similar to ERFNd for most latitudes, except in the Southern Ocean where there was a much larger ERFLWPic contribution. ERFfc contributions were mostly between around 50 % and 75 % of the ERFNd contribution, but again in the Southern Ocean there was a larger ERFfc contribution than ERFNd contribution. This was also true in the polar regions. Globally, the overall contribution from adjustments (ERFfc and ERFLWPic) was 0.92 W m−2, compared to an ERFNd contribution of 0.52 W m−2, so that in this model more aerosol forcing is coming from the more complicated adjustment processes, highlighting the importance of evaluating these processes in more detail. Here, we perform a similar analysis but with a different model. This will allow the two models to be compared in terms of how they produce their aerosol forcing. A very different breakdown between the models would mean that one of the models was incorrect and allow the basis for more detailed evaluation. We also focus on the NA region and on surface aerosol forcing rather than top-of-the-atmosphere (TOA) forcing due to the potential importance of aerosol forcing for the climate variability via sea surface temperature changes there. The focus on the NA also allows a more detailed look at the processes than is possible from a global study. We also go beyond the study of Mülmenstädt et al. (2019) by also characterizing the cloud regimes and the changes in cloud regimes for which aerosol forcing predominantly occurs according to the model; this will then allow (in future work) these regimes to be comprehensively evaluated against observations and also against high-resolution models. High-resolution modelling needs to be targeted to a smaller selection of regimes given the high computational cost. In addition, observations will be used to evaluate the general model cloud properties (spatial positioning, areal coverage, thickness, droplet concentrations, etc.) since these will affect the resulting aerosol forcing.

2 Methods

2.1 Definition of liquid water path

Throughout this paper, the term “LWP” refers to the all-sky value (i.e. including both the cloudy and clear-sky portions of model grid boxes or observed regions) and LWPic refers to the in-cloud liquid water path, which is that from the cloudy regions only. It is assumed here that LWP =fcLWPic (e.g. as also in Seethala and Horváth2010).

2.2 Model details

We examine the aerosol–cloud interactions in the UKESM1-A model, which is the atmosphere-only version of the coupled UKESM1 (UK Earth System Model), which was submitted as part of the sixth Coupled Model Intercomparison Project (CMIP6). UKESM1 is based on the HadGEM3-GC3.1 physical climate model (Williams et al.2017; Kuhlbrodt et al.2018) but in addition couples several Earth system processes (Sellar et al.2019). These additional components include the MEDUSA ocean biogeochemistry model (Yool et al.2013), the TRIFFID dynamic vegetation model (Cox2001) and the stratospheric–tropospheric version of the United Kingdom Chemistry and Aerosol (UKCA) chemistry model (Archibald et al.2020). This version of the UKCA allows a more complete description of atmospheric chemistry compared to HadGEM3-GC3.1. For example, the latter uses an offline climatology for oxidants, whereas in the UKESM1 oxidants are treated explicitly.

The UKESM1-A atmosphere-only version differs from the UKESM1 in that it does not include the ocean and sea-ice models, MEDUSA and TRIFFID. Instead, the UKESM1-A configuration uses observed sea surface temperatures and sea-ice concentrations provided by the Program for Climate Model Diagnosis and Intercomparison (Taylor et al.2000, last access: 11 December 2020). Vegetation (vegetation fractions, leaf area index, canopy height) and surface ocean biology fields (dimethyl sulfide and chlorophyll ocean concentrations) are prescribed from a member of the UKESM1 CMIP6 historical ensemble (Sellar et al.2019). This ensures that the prescribed vegetation and ocean biological fields mirror those simulated by the TRIFFID dynamic vegetation and MEDUSA scheme in the coupled historical run. The vegetation fractions are prescribed from time-varying annual means; leaf area index, dimethyl sulfide (DMS) and chlorophyll concentrations are monthly values from the time means of the 1979–2014 period; and the canopy heights are an overall time mean for 1979–2014. All other emissions of gases and particles from sea and land surfaces are prescribed from the CMIP6 inventories; see Mulcahy et al. (2020) for details of the specific implementation for the UKESM1.

The atmospheric component is the GA7.1 atmospheric configuration of the Unified Model (UM). Full details of this can be found in Walters et al. (2019) and Mulcahy et al. (2020). Here, we only describe in detail the features that are more relevant to aerosol–cloud interactions. We use an N96 horizontal resolution, which is 1.875 × 1.25 (208 × 139 km) at the Equator. A total of 85 vertical levels are used between the surface and 85 km altitude with a stretched grid such that the vertical resolution is 13 m near the surface and around 150–200 m at the top of the boundary layer. We chose this resolution since it is the same as that used for long climate runs in the CMIP6 model intercomparison and is the same horizontal resolution as that used in the B12 study.

Aerosol number concentrations are treated prognostically with the Global Model of Aerosol Processes (GLOMAP) multi-modal scheme (Mann et al.2010, 2012), which uses five log-normal aerosol size modes and includes sulfate, sea salt, black carbon and organic matter chemical components. These aerosol components are treated as being internally mixed within each size mode. Mann et al. (2010) and Mann et al. (2012) provide further details with some small changes for the implementation within UKESM1-A described in Mulcahy et al. (2020). Mineral dust is simulated separately using the CLASSIC dust scheme (Woodward2001; Mulcahy et al.2020).

Shallow, middle and deep convection are parameterized separately to other cloud types (see Walters et al.2019, for details). The parameterizations do not take into account aerosol, nor droplet concentrations, and they use their own simplified microphysics scheme. Therefore, the representation of aerosol forcing is incomplete. For clouds that are not shallow, middle or deep convection (termed large-scale clouds), the UKCA-Activate scheme is used to calculate cloud droplet concentrations from the aerosol size distribution using the West et al. (2014) scheme based on the parameterization of droplet activation in Abdul-Razzak and Ghan (2000). Cloud droplet concentrations at cloud base are replicated vertically throughout contiguous columns of clouds. The cloud droplet activation scheme is a diagnostic scheme since it is run on each model time step without consideration of how many cloud droplets were present before. The cloud droplet number concentration is then passed to the radiation and microphysics schemes.

The large-scale cloud microphysics is a single-moment scheme, such that the mass of liquid water, but not the cloud droplet number, is advected by the model and retained in memory between model time steps. It is based on Wilson and Ballard (1999) but with improvements to the warm rain parameterizations suggested by Boutle et al. (2014), which include bug fixes, a treatment of rain fraction that is consistent with the prognostic rain formulation, a switch to the Khairoutdinov and Kogan (2000) parameterization for autoconversion and accretion, and a bias correction for the latter processes to deal with subgrid variability of clouds and rainwater. Relative to Wilson and Ballard (1999), there is also an improved treatment of drizzle rates (Abel and Shipway2007) and a prognostic treatment of rain that allows the three-dimensional advection of precipitation. The introduction of the latter required modifications to be made to the aerosol wet scavenging processes as described in Mulcahy et al. (2020). The bulk properties (cloud fraction, cloud liquid water content, vertical overlap, etc.) of large-scale clouds are parameterized using the prognostic cloud fraction and prognostic condensate (PC2) scheme (Wilson et al.2008a, b) with modifications described in Morcrette (2012). The atmospheric boundary layer is parameterized using the turbulence closure scheme of Lock et al. (2000) with modifications described in Lock (2001) and Brown et al. (2008).

There are some differences in the treatment of aerosols in UKESM-A and the physical climate model (HadGEM-GC3.1) primarily related to natural aerosols, aerosol chemistry and the prescription of anthropogenic SO2 emissions (see Mulcahy et al.2020, for details).

2.3 Simulation details

The model is run from 1 March 2009 to 28 March 2010. The first 27 d were used to spin up the aerosol and chemistry fields, leaving 1 year of data for analysis. The time period was chosen to allow comparisons with relevant field campaigns that will be performed in future work. Two parallel simulations were performed: one using pre-industrial (PI) CMIP6 aerosol emissions from the year 1850 and the other using present-day (PD) emissions (i.e. the CMIP6 emissions corresponding to the simulation period). Both simulations were nudged every 6 h to ERA-Interim horizontal wind fields between  2277 and  47 251 m (applied between the 18th and 76th model levels from the surface). The nudging ensures that the meteorology in the two runs is very similar, thereby allowing cloud radiative effects due solely to the aerosol perturbation to be calculated. Following the recommendations of Zhang et al. (2014), we do not nudge the temperature field in order to allow fast-acting local responses to aerosol-induced temperature changes (such as those from precipitation suppression, ARI and semi-direct radiative effects). Instantaneous model fields are output every 27 h, which allows the sampling of the diurnal cycle and more complex output analysis than is possible using time-averaged data.

The results in the main body of the paper are based on meteorology and emissions from the period of 28 March 2009 to 28 March 2010. It is possible that the results presented vary depending on the chosen year since meteorology, cloud fields, etc. vary from year to year. To address this, we have also run the PI and PD simulations for an additional year for the period of 28 March 2010 to 28 March 2011. We have compared selected key results from Sect. 3.2, which are shown in Appendix G. Very similar results are found using the alternative year, which demonstrates that our results are robust and not sensitive to the chosen year of meteorology.

2.4 Surface forcing calculations and forcing partitioning

In this paper, we are concerned with the shortwave aerosol ERF at the surface since we are interested in the effects on SSTs and NA climate variability. To separate the aerosol ERF into ARI and ACI components, we use output from the triple calls to the radiation scheme that are made by the model every time step. One call calculates the surface SW fluxes taking into account all radiatively active components including aerosols and clouds (designated here as SWaerosol + cloudy), one call calculates the fluxes with the background aerosol removed (i.e. the aerosol outside of clouds and interstitial aerosol; aerosol Twomey and adjustment effects on the clouds are still included here; SWclean+cloudy), and one call calculates the fluxes with both the background aerosol and clouds removed (SWclean+clear). These are detailed in Table 1. The instantaneous changes in SW due to ARIs and ACIs for a given grid box are then calculated as

(1) Δ SW ARI = SW aerosol + cloudy - SW clean + cloudy Δ SW ACI = SW clean + cloudy - SW clean + clear .

The corresponding instantaneous radiative forcings due to anthropogenic aerosols (ERFARI and ERFACI) can then be calculated from simulations with PI and PD aerosol emissions:


Table 1Details about the aerosol and clouds that are used in the three calls to the radiation scheme performed on each time step.

Download Print Version | Download XLSX

This is similar to the technique used in Ghan (2013) and Jiang et al. (2016). We also decompose the surface ACI aerosol forcing into contributions from changes in Nd, LWPic and fc between the PI and PD. This is achieved using an offline calculation of the surface SW fluxes, which allows the magnitude of each term to be assessed individually. The approach is described in Grosvenor et al. (2017) and in Appendix F. We also note that TOA fluxes could also be decomposed using this technique, but we focus on the surface forcing here.

3 Results

3.1 Model evaluation against satellite observations

Here, we evaluate the PD simulation against satellite observations but limit our analysis to ocean regions due to the lesser reliability of satellite products over land. The motivation for the model evaluation is that without a good representation of cloud properties and spatial distribution it is likely that the model aerosol forcing will be incorrect. For example, if the simulated clouds have a cloud fraction that is too high, then a cloud albedo perturbation due to aerosol would be larger than it is in reality. Placement of clouds in the wrong locations could affect forcing via the differing incoming solar insolation, or it could affect their chances of interacting with sources of aerosol. Where possible, we use the COSP (Cloud Feedback Model Intercomparison Project Observation Simulator Package; Bodas-Salcedo et al.2011) satellite simulator model output for the relevant satellite in order to get a fairer comparison between the model and satellite. This is particularly important when comparing cloud fractions since the cloud detection threshold (e.g. in terms of optical depth) varies between the different satellite instruments and needs to be consistent with the threshold LWC used to define a cloud in the model. Furthermore, COSP accounts for the effect of overlying layers of clouds, which affects the observation of underlying clouds.

Following Gustafson and Yu (2012), model biases are quoted in terms of the normalized mean bias factor (NMBF) and the normalized mean absolute error factor (NMAEF), which are both expressed as a percentage, in order to provide unbiased metrics that are symmetric about zero. These quantities are similar to the mean percentage bias and the root mean square error (RMSE) except that account is taken of whether the model is greater or lower than the observations. In addition, the NMAEF does not involve squaring the error, which can lead to an exaggerated influence from larger biases. Appendix B provides the definitions for these metrics. Table 2 lists these bias metrics (along with the spatial correlation coefficient r, between the model and the observations) for the different regions considered (see next section) and for the different cloud variables. It should also be borne in mind that, as well as being due to issues with the representation of clouds in the model, differences in cloud properties between the model and satellite might also be due to cloud adjustments to aerosol that are too strong or weak, and that it is difficult to distinguish between the two.

3.1.1 Cloud fraction evaluation

The distribution of time-mean low-altitude fc from the model is compared to CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observation) satellite lidar observations in Fig. 1. Version 3 of the CALIPSO-GOCCP (GCM Oriented Cloud Calipso Product; Chepfer et al.2010) is used, which is available at, last access: 11 December 2020. The model has a good spatial correlation with the satellite (r=0.83) over the region shown, with both model and observations showing higher cloud fractions in the northern part of the NA where the time-mean values can reach around 70 % especially towards the west near the coast of Canada and north of Iceland and Scandinavia. More than 40 %–50 % of the low-altitude clouds in these regions are stratocumulus (Wood2012). The model has a positive bias of around 30 %–40 % in the region off the eastern Canadian coast and along the east coast of the US. There is also a band of higher cloud fractions (around 20 %–30 %) down the west coast of Africa in both the model and observations; this region is also associated with a high occurrence of stratocumulus (40 %–60 % of the low clouds according to Wood2012). Model biases are minimal in these regions. The model and observations both show lower cloud fractions in the southwest region of the domain.

The model captures the main pattern of low-altitude cloud fraction well, giving us some confidence that it provides a realistic representation of the broad cloud types and locations in the NA. However, the model uses observed SSTs, and the winds above the boundary layer are nudged to reanalyses, which will help the model to replicate the correct cloud patterns to some degree. But even with the correct SSTs and wind fields, accurate simulation of the spatial pattern of fc is a strong test of the model boundary layer and cloud schemes, in particular the correct vertical thermodynamical structure and entrainment, the corresponding cloud fraction response, etc. Furthermore, the free-running ocean-coupled version of the model (UKESM1; Sellar et al.2019), which has no wind nudging and predicts the SSTs itself, also reproduces the cloud pattern and amount well (Robson et al., 2020). This also implies that the coupled model and the nudged model used here exhibit similar cloud regimes and gives more confidence that the results in this paper apply to coupled models.

We define a subregion (marked in Fig. 1 and denoted northern North Atlantic, or northern NA) to represent the region where the stratocumulus cloud fraction is large. We also define a second region immediately to the south of the northern NA region (southern North Atlantic, or southern NA; see Fig. 1) where the annual mean stratocumulus cloud cover is lower (Wood2012), indicating more broken stratocumulus and/or a lower frequency of occurrence of stratocumulus. The overall low-altitude fc area-weighted biases for these regions were 5.1 % and 23.9 % for the northern and southern NA regions, respectively (see Table 2).

Figure 1Low-altitude (CTP  680 hPa) cloud fraction evaluation. (a) Model; (b) CALIPSO satellite data; (c) model bias. For the model, we are using the COSP simulator for CALIPSO to deal with the satellite cloud detectability threshold and overlying layers.

Figure 2As for Fig. 1 except for mid-altitude clouds (440  CTP < 680 hPa).

For mid-altitude clouds (Fig. 2), the spatial pattern is again captured well by the model (r=0.80). Percentage biases are generally low in the northern NA region with a mean model bias of 11.0 %. There are larger negative biases in the southern part of the domain with a mean bias of 25.8 % for the southern NA region; however, the observed cloud fractions are very low in those regions, making higher percentage biases more likely. There are also positive biases of up to around 30 % in the stratocumulus region north of Scandinavia.

The modelled high-altitude cloud fraction is generally biased high over the ocean, with mean biases of 7.3 % for the northern NA and 4.5 % for the southern NA, although the spatial pattern is again good (r=0.79 for the whole region). Positive biases are particularly evident in the region north of around 60 N, where biases of up to 90 % occur. Since low-altitude clouds are likely to be the most important for aerosol–cloud interactions, we will leave an investigation of the biases in the mid- and high-altitude clouds to future research.

Figure 3As for Fig. 1 except for high-altitude clouds (CTP < 440 hPa).

Table 2Model evaluation statistics for the various subregions using time-averaged data. The normalized mean bias factor (NMBF) and the normalized mean absolute (NMAEF) error factor statistics are used following Gustafson and Yu (2012); see Appendix B for definitions. r is the spatial correlation coefficient between the model and observed time averages. All values are area weighted to account for the variation in area of the model grid boxes.

Download Print Version | Download XLSX

3.1.2 LWP evaluation

Liquid water path is important for cloud optical depth and hence cloud albedo. For example, for a fully overcast cloud with a liquid water content that increases with height in a manner consistent with adiabatic uplift, the optical depth (τc) is proportional to LWPic5/6Nd1/3. Therefore, a given relative change in LWPic produces more than twice the relative change in τc that the same relative change in Nd does. LWPic is also important in terms of cloud physics since it helps to determine rain rates and droplet sizes.

Microwave satellite instruments such as AMSR-E (Advanced Microwave Scanning Radiometer for Earth Observing System; e.g. Wentz and Meissner2000), whilst only being able to retrieve over ocean surfaces, probably represent the most accurate of the available retrievals of LWP since they are not subject to the biases associated with retrievals of LWP that use a combination of visible and shortwave-infrared wavelengths (e.g. Moderate Resolution Imaging Spectroradiometer, MODIS; Salomonson et al.1998), they give a better representation of the total column LWP than from such retrievals, they are not affected by the presence of ice, and they are available during both the daytime and nighttime (O'Dell et al.2008; Elsaesser et al.2017). Note that LWP from AMSR-E is the all-sky LWP and so includes the zero values from clear regions and the LWP contributions from cloudy regions (as does the LWP output from the model), whereas MODIS retrieves the in-cloud LWP, LWPic. However, AMSR-E retrievals are still subject to biases; Lebsock and Su (2014) suggest that the main cause of bias is the presence of rain (and the inability to directly detect whether rain is present). Therefore, as suggested in Elsaesser et al. (2017), we place more confidence in the microwave LWP observations when the ratio of the all-sky LWP to total liquid water path (TLWP) (the sum of LWP and RWP, where RWP is rainwater path) is high; we denote the ratio as fLWP. A caveat here, though, is that the partitioning of LWP and RWP from microwave instruments, and hence the estimate of fLWP, is itself uncertain.

Figure 4Time average of fLWP (ratio of LWP to the sum of LWP and RWP). For the model, LWP and RWP contributions from both the large-scale cloud scheme and the convective parameterization are included. The RWP for AMSR-E is calculated from the retrieved LWP and rain rate by inverting the retrieval algorithm, as described in Elsaesser et al. (2017). The same raindrop size distribution and fall speed relationship as for the model is used. Both daytime and nighttime AMSR-E overpasses are used.

Figure 5Time-mean in-cloud LWP (LWPic) model evaluation for both day and night overpasses. The bias plot on the right includes a contour of the 0.9 value of fLWP; see Fig. 4 for the full map of fLWP for reference.

Given the uncertainties when rain is present, we also initially only use the model LWP and not the RWP. Furthermore, we only consider the model LWP from the large-scale cloud scheme and not that from the convective parameterization. Convective LWP will be mostly associated with heavily raining environments and so the observations in such regions are again more uncertain. However, we examine and discuss how RWP and convective LWP and RWP change the model evaluation in Appendix D. The AMSR-E instrument is aboard the Aqua satellite, which has local overhead overpass times of 01:30 and 13:30 LT, and we use data from both overpasses. Therefore, when comparing to AMSR-E for LWP and RWP, the model local times within 3 h on either side of both of these times are used. This is important for LWP and RWP because cloud water content can have a large diurnal cycle.

Figure 4 shows fLWP from the model and the AMSR-E satellite along with the model bias relative to the retrieval. fLWP is quite high throughout the region in both the model and the observations, with values over the ocean generally larger than around 0.8. The model bias is generally small (< 10 %), except off the east coast of the US where there is an overestimate, although with biases mostly less than 0.1 ( 15 %). This region corresponds to a region of high LWP (see Fig. 5) and might suggest that too much rain is produced by the model, although, given that fLWP is low, it is also a region where the satellite estimate of fLWP is likely to be more uncertain. The largely good agreement between the model and satellite for fLWP provides some confidence in the ability of the model to represent cloud and rain formation processes, but insofar as we assume the model to be realistic, it also lends confidence to the satellite estimates of this quantity for this region, which, as explained above, is also uncertain.

Figure 6As for Fig. 5 except both the model and satellite data have been filtered before time averaging to only include data points for which fLWP is greater than 0.99. This quantity is denoted as LWPic0.99.

LWPic values were estimated by dividing the time-mean all-sky LWP by the time-mean low-altitude cloud fraction from CALIPSO. The same was done for the model using the time-mean CALIPSO low cloud fraction from the COSP satellite simulator. For reference, the evaluation of the grid-box mean LWP (including cloudy and clear regions) is shown in Appendix D.

Figure 5 evaluates the time-mean LWPic with no filtering using fLWP. The bias plot in Fig. 5 shows the 0.9 value of the time-mean satellite fLWP contour to highlight regions where the satellite data are likely to be more reliable. The spatial patterns of the model and satellite have a spatial correlation r value of 0.68. There is a negative model bias off the east coast of Florida, although there is generally a lot of rain present in this region, as indicated by the fLWP contour and Fig. 4, such that the satellite measurements are more uncertain. There is also a negative model bias between Greenland and Canada and on the east side of Greenland where the fLWP is high. The grid-box mean LWP is also biased low there (Fig. D1), indicating that the model clouds are too thin (rather than the cloud fraction too low; see also Fig. 1). Positive biases occur off the west coast of Africa at around 18 N; however, the LWPic is very low in this region. Further positive biases occur in the stratocumulus to the north of Scandinavia. Again, these are not associated with cloud fraction biases (see Figs. D1 and 1), suggesting model stratocumulus clouds that are too thick.

Figure 6 compares LWPic after filtering both the model and satellite data (before time averaging) to only include data points for which fLWP>0.99 (hereafter LWPic0.99). This selects cloud scenes that are likely to be non-precipitating and non-convective clouds, which are likely to be either stratocumulus or shallow cumulus clouds. The regions of stratocumulus north of Scandinavia and around Iceland where there were positive biases with no filtering for fLWP>0.99 now have even larger LWPic biases, again indicating that the model low clouds may to be too thick. The negative biases west of Greenland remain. Since there is generally little RWP in this region and the fLWP filtering has been performed, this indicates that the LWPic evaluation in this region is likely to be robust. There are now negative biases extending from the coast of Florida up to near the UK, suggesting that clouds there are too thin, although instrumental uncertainties may still be high here given the uncertainty in the satellite estimate of fLWP and the larger amount of RWP that occurs here.

The LWPic values from the AMSR-E observations in Fig. 6 show remarkably little spatial variability, suggesting that non-precipitating clouds tend to be fairly uniform across broad regions with an LWPic value of around 60–100 g m−2. Cloud fraction (Fig. 1) and Nd (Fig. 7) vary quite widely across the same region, potentially indicating a physical process that maintains a constant LWPic despite a varying aerosol environment and cloud regime.

Figure 7Time-mean model Nd evaluation vs. the MODIS satellite. The model and MODIS data are restricted to data points for which the grid-box mean cloud top height is ≤3.2 km and for which the liquid cloud fraction is ≥80 %. The COSP liquid cloud fraction is used for the model screening. The arrangement of the panels is as for Fig. 1.

3.1.3 Cloud droplet concentration evaluation

Cloud droplet concentration (Nd) is an important quantity because it generally represents the first step in the chain of processes by which aerosols affect clouds. Nd gives some indication of the number of cloud condensation nuclei (CCN) that were available to produce clouds as well as other factors such as updraft speed, droplet collision coalescence, droplet scavenging by rain, cloud evaporation, etc. Thus, an evaluation of model Nd can give an idea of how well the model captures a range of processes.

Here, we evaluate Nd using a 1×1o resolution data set calculated from MODIS retrievals of τc and re. Two-dimensional fields of Nd are derived by the retrieval since it is assumed that Nd is constant throughout the depth of the cloud, which has been shown to be a good approximation by aircraft observations of stratocumulus (Painemal and Zuidema2011). Details of the retrieval are given in Appendix A. For the model, two-dimensional Nd fields were obtained from the instantaneous 3D Nd fields by calculating a weighted vertical mean Nd, with the LWC on each level used for the weights. This ensures that the levels with the most LWC contribute most to the average Nd, which is similar to what is assumed in the MODIS retrieval since most of the re signal comes from near cloud top where the LWC is assumed to be the largest and the Nd calculation is a strong function of re. It also avoids contributions from very thin clouds that would not be detected by MODIS. Only data points for which the liquid cloud fraction is larger than 80 % and for which the mean cloud top height is below 3.2 km were included for the satellite Nd calculation in order to help exclude satellite retrieval errors (see Grosvenor et al.2018b). The same filtering was applied to the model to minimize sampling errors; the COSP MODIS liquid cloud fraction and a calculation of model cloud top height were used for this. The satellite data set was regridded to the model grid before comparison.

Figure 7 shows that the modelled and observed Nd values are largest near the continents and that Nd decreases towards the middle of the Atlantic Ocean. This fits with the idea that CCN are scavenged during eastward transport (Wood et al.2017). There are also likely to be some dilution effects as distance from the sources increases. Other factors may also influence the spatial Nd pattern such as changes in boundary layer height, predominant cloud type, other meteorological factors, etc. CCN concentrations also increase close to the European and African source regions, consistent with the prevailing northerly wind, which transports European and African pollution down the coast. Overall, the model produces a reasonable spatial pattern with a spatial correlation coefficient of 0.67 for the whole region. However, the model has a tendency to overestimate Nd over the southern part of the North Atlantic region, off the east coast of the UK and over Europe. Conversely, the model underestimates in the northern part of the domain. The NMBF is 20.6 % for the northern NA region where stratocumulus dominates compared to 21.5 % for the southern NA region where the clouds are more broken.

Figure 8Time-mean model SWTOA evaluation vs. the CERES satellite. The monthly mean CERES-EBAF data product is used here; this product uses data from both the Terra and Aqua satellites as well as geostationary satellites in order to approximate averaging across the diurnal cycle. The arrangement of the panels is as for Fig. 1.


3.1.4 Shortwave top-of-the-atmosphere flux evaluation

Clouds and the Earth's Radiant Energy System – Energy Balanced and Filled (CERES-EBAF) data were taken from, last access: 11 December 2020 and are the monthly averaged product of observed TOA for which the TOA net flux is constrained to the ocean heat storage.

Figure 8 shows the evaluation of the time-mean top-of-the-atmosphere SW radiative fluxes (SWTOA) vs. the CERES satellite. The spatial pattern of the model matches the observed pattern well with a correlation coefficient of 0.89. The model biases are mostly small, with positive biases that are <20 % over all of the oceanic regions except south of the Equator where positive biases of up to <35 % occur. The NMBF bias for the whole region (including land) is 2.4 % and it is 5.8 % and 1.7 % for the northern NA and southern NA regions, respectively.

Figure 9Estimates of the time-mean SWTOA flux (a, b; W m−2); (a) calculated using the time-mean modelled cloud properties (for comparison with Fig. 8a); (b) calculated using the time-mean observed cloud properties (for comparison with Fig. 8b). Panel (c) shows the sum of the estimated contributions from the individual cloud property biases (see Fig. 10) expressed as a percentage of the observed CERES SWTOA flux (for comparison with Fig. 8c).

Figure 10Estimated contributions from the model biases in Nd (a; note the smaller colour bar range), LWPic (b) and fc (c) to the SWTOA model bias (Fig. 8c). Values are expressed as a percentage of the observed CERES SWTOA field for consistency with Fig. 8c. Figure 9c shows the sum of these three contributions also as a percentage of the observed CERES SWTOA.

The SWTOA biases can be caused by a combination of biases in fc, TLWPic and Nd. To estimate the relative contributions of these different biases individually to the SWTOA bias, we used a similar calculation to that described in Sect. 2.4 but applied to TOA fluxes, following the technique described in Grosvenor et al. (2017). Firstly, the modelled time-mean SWTOA field was reconstructed by using the modelled time-mean values of the COSP CALIPSO fc, LWPic (with no filtering using fLWP) and Nd (filtered as in Sect. 3.1.3) as inputs into the SWTOA offline calculation. This was also done using the time-mean observed fields from CALIPSO, AMSR-E and MODIS. The results are shown in Fig. 9 and can be compared to Fig. 8a and b. Comparisons are not made over land since the higher land albedo was not taken into account in the calculations. In both cases, the spatial pattern of the calculated SWTOA (Fig. 9a and b) over the ocean correlates well with that of the actual SWTOA (Fig. 8a and b).

The SWTOA values from the calculations are somewhat smaller than those actually modelled or observed. If considering only oceanic regions, the underestimates for the model off the east coast of Canada in the NA stratocumulus region are perhaps of greatest concern, since this is where the SWTOA values and their biases were highest. The underestimates are likely due to the approximations and assumptions that have been made with this approach such as using the time-mean shortwave flux and cloud fields. However, the agreement is sufficient for the purposes of estimating the relative contributions from biases in individual cloud properties to the SWTOA bias.

Figure 11The percentage of the sum of the absolute biases contributed by the biases in Nd (a), LWPic (b) and fc (c). See Eq. (3). Note that the absolute values of the percentages for each grid box add up to 100 %.

Next, perturbations to SWTOA were calculated by applying the model biases in the individual cloud fields (fc, LWPic and Nd) to the observed cloud values on a one-at-a-time basis (Fig. 10). They are expressed as a percentage of the observed time-mean SWTOA field. The sum of these perturbations, again expressed as a percentage, is shown in Fig. 9c and is intended for comparison to the actual model SWTOA bias shown in Fig. 8c. The spatial pattern of the combined bias estimate matches the actual bias well, although there are some regions of negative bias that are not present in the actual bias. These are in regions of low SWTOA and low fc and so are more prone to error and likely less important for the overall SWTOA. In the regions of positive model SWTOA bias, the estimate is a little lower than the true model bias, but again the agreement should be sufficient to compare the relative contributions from the individual cloud field biases.

From Fig. 10, it is clear that perturbations of the magnitude of the fc model biases have a very large impact on the SWTOA and are the main contributor to model SWTOA biases in most regions. Large negative contributions occur to the SWTOA bias in the southern part of the domain due to the negative cloud biases there (Fig. 1). Smaller positive contributions from Nd and LWPic biases also occur in this region to give less overall estimated SWTOA bias (Fig. 9c) in agreement with the small SWTOA model bias (i.e. directly from the model output; Fig. 8c). This suggests that some cancellation of biases in fc, LWPic and Nd is occurring here resulting in the observed low SWTOA bias. Large positive contributions to the SWTOA bias from fc biases occur in the stratocumulus region in the western part of the northern NA. These are also partially cancelled by negative Nd and LWPic biases, but the overall SWTOA is positive in this region (Fig. 8c). LWPic (note the different colour scale) has only a minor effect north of Scandinavia (contributing to positive biases) and west and east of Greenland (negative contributions). The Nd biases contribute even less to the SWTOA model bias than LWPic biases; positive contributions from Nd biases are seen south of 36 N and negative contributions north of there.

Next, an attempt to quantify the relative percentage contributions (Px) of the individual cloud field properties (denoted as x, where x is either fc, LWPic or Nd) to the total SWTOA bias is made using the following equation:

(3) P x = 100 × Δ SWTOA x | Δ SWTOA f c | + | Δ SWTOA LWP ic | + | Δ SW N d | ,

where ΔSWTOAx is the perturbation in SWTOA due to the bias in cloud field property x. Figure 11 shows that fc biases dominate over nearly all regions south of around 54 N. However, there are some regions where LWPic biases dominate; there are large positive relative contributions from LWPic in the region north of Scandinavia and some large negative contributions surrounding Greenland. Nd contributions are important in the diagonal band stretching from the southwest of Iceland up to Svalbard.

Overall, this analysis shows that positive biases in fc are important in explaining the positive SWTOA biases east of the US and Canada, while positive LWPic biases are important off the northern coast of Scandinavia.

3.2 Aerosol forcing

We now examine maps of the temporal mean aerosol forcings following Eq. (2). Note that the calculations were performed using the instantaneous fields from the 27-hourly output, which ensures that the diurnal cycle is sampled evenly throughout the course of the simulation.

Figure 12Time-mean ERFARI (a) and ERFACI (b) forcing. The maps have been averaged over 3×3 grid boxes to help reduce noise.

Figure 12 shows maps of ERFARI and ERFACI. The magnitude of ERFACI is generally larger than that of ERFARI for oceanic regions. For example, in the northern NA box, the mean ERFARI is 0.50 W m−2 and ERFACI is 3.1 W m−2 (see Table 3). For the southern NA box, ERFARI and ERFACI are 1.0 and 1.8 W m−2, respectively. The dominance of ERFACI over the NA region is in agreement with B12 (Fig. 4b of that paper). The overall spatial pattern of ERFACI is also similar to that in B12 with a band of large negative forcing running across the Atlantic from west to east between around 25 and 50 N, and another region of negative forcing following the west coast of Africa. Some differences are apparent, though; for example, the negative forcing is smaller in magnitude southwest of the UK in our simulations and larger in magnitude near Africa south of the Equator. However, over the whole NA domain, ERFARI is larger in magnitude than ERFACI (1.9 vs. 1.7 W m−2) due to the fact that ERFARI is usually larger over land. The dominance of ERFARI in the southern regions of the domain will also have more influence on the area-weighted average.

Table 3Temporal and spatial means of various quantities for the whole NA domain and the various subregions. Values are shown for ERFARI and ERFACI, and the percentage changes in Nd, LWPic and fc between the PI and PD runs. Percentage fc changes are also shown for the runs with no aerosol effect on the rain autoconversion process. The other columns show results from the offline estimates of the ACI forcing. ERFACI,all is the estimated forcing due to simultaneous PI to PD perturbations of Nd, LWPic and fc. ERFNd, ERFLWPic and ERFfc are the forcing estimates from one-at-a-time perturbations of Nd, LWPic and fc, respectively, and ERFACI,sum is the sum of these three forcing estimates. ERFACI,macro is the forcing contribution from the cloud macrophysical changes, i.e. the sum of ERFLWPic and ERFfc. All values are area weighted to account for the variation in area of the model grid boxes.

Download Print Version | Download XLSX

3.2.1 Changes in cloud properties from PI to PD and their contributions to the ACI forcing

Figure 13 shows the time-mean changes in cloud properties (Nd, LWPic and fc) between the PI and PD simulations (PD minus PI as a percentage of PI). Considering oceanic regions, percentage increases in Nd are greatest off the east coast of the US and Canada (to the south of Newfoundland), off the SW coasts of Spain and Portugal and west Africa, and in the region north of Scandinavia. There are increases across the entire Atlantic running from the US to west Africa, but the magnitude decreases moving east (except close to the European–African west coast) likely reflecting the decreasing influence of pollution from the east coast of the US. The spatial pattern of Nd change matches the spatial pattern of the change in column sulfate aerosol mass fairly well (Fig. C3) suggesting that changes in sulfate aerosol are the main cause of the Nd changes. The exception is the region stretching from southern Portugal down the west coast of north Africa and also across the Atlantic south of around 20 N, where the Nd changes coincide with changes in black carbon (BC) and organic matter (OM) column mass.

LWPic changes are generally quite noisy, but the largest changes over the ocean occur west of Greenland and in a diagonal band across the Atlantic in a similar way to the Nd changes except further to the north. The changes in fc are largest in the southern regions of the domain below around 35 S. The responses of LWPic and fc can be termed macrophysical responses (or also adjustments), since they affect bulk cloud properties, as opposed to the response of Nd that is termed a microphysical response. The macrophysical responses are mostly due to the suppression effects of aerosols on rain rates via the autoconversion process. This is demonstrated in Figs. E1, E2 and Table 3, which show the impact of preventing aerosol from affecting the autoconversion process (see Appendix E for details). We hypothesize that the changes in LWPic over the northern NA regions and lesser changes in fc reflect the presence of stratocumulus clouds with very large areal coverage in the north that can only respond to the suppression of rain by thickening (i.e. more LWPic) rather than by increasing fc. In the southern NA, the cloud coverage is much lower due to the presence of broken stratocumulus and cumulus clouds, which respond in a different way to the rain suppression, namely by increasing their coverage (fc) more than their thickness (LWPic). We note that in reality clouds can also respond to enhanced aerosol by increasing cloud top entrainment, which can reduce LWPic and fc (Ackerman et al.2004; Bretherton et al.2007; Hill et al.2009). However, this mechanism is currently not included in the model.

Based on the maps of the cloud property changes alone, it is difficult to judge the relative contributions of each type of change to the forcing. Thus, we now make a quantitative estimate of this using offline radiative calculations as described in Sect. 2.4 and Appendix F. Figure 14a shows the ACI forcing estimated using the offline method following Eq. (F13) for the case where all of the cloud variables have been perturbed from their PI values to their PD values (ERFACI,all). Figure 14b shows the sum of the contributions calculated by separately perturbing the individual cloud properties in one-at-a-time tests (ERFACI,sum). These can be compared to the forcing diagnosed from the full model outputs as generated by the online radiation code (Fig. 12b). The general pattern and magnitude of the actual and estimated forcings are very similar. The fact that ERFACI,all and ERFACI,sum are similar indicates linearity in the effect of the perturbation of the individual cloud properties, such that they are almost directly additive. The mean forcing over the northern subregion (northern NA) is 3.2 W m−2 for ERFACI,sum, 3.0 W m−2 for ERFACI,all and 3.1 W m−2 for the actual forcing. In the southern subregion (southern NA), the mean forcings are 1.9, 2.0 and 1.8 W m−2 for ERFACI,sum, ERFACI,all and the actual forcing, respectively. Overall, the match is very good, suggesting that the estimation technique is sound and that it is likely to be useful for estimating the contributions from the different cloud property changes to the forcing.

Figure 15 shows the contribution to the surface ACI forcing from the changes in Nd, LWPic and fc following Eq. (F15). Percentage contributions (Px) from the individual cloud fields (denoted as x, where x is either Nd, fc or LWPic) to the sum of the absolute contributions are calculated in a similar way to Eq. (3):

(4) P x = 100 × ERF x | ERF N d | + | ERF LWP ic | + | ERF f c | ,

where ERFx is the ACI forcing contribution due to the change in cloud field x. These values are mostly negative since the ACI forcing is negative, but positive contributions and thus positive ERFx are possible. Figure 16 shows the results.

From the summary of the area means in Table 3, we see that for the northern NA region the area-mean contributions from changes in Nd and LWPic (ERFNd and ERFLWPic) are similar with values of 1.4 and 1.2 W m−2, respectively. These are larger than the fc contribution (ERFfc) of 0.54 W m−2. In contrast, for the southern NA region ERFfc is 1.2 W m−2, which is considerably larger than the Nd or LWPic contributions (0.63 and 0.23 W m−2, respectively). Overall, for the whole NA region, ERFfc is 0.74 W m−2, which is similar in magnitude to that from the Nd changes (0.67 W m−2) with ERFLWPic a little smaller at 0.41 W m−2.

The sum of the macrophysical contributions from fc and LWPic (termed ERFACI,macro) is the dominant contribution to aerosol forcing in the NA region providing 63.9 % of ERFACI,sum (= ERFNd+ ERFLWPic+ ERFfc). This is important because these changes are likely to be associated with a higher model uncertainty than pure cloud albedo effects (due to Nd changes alone). The macrophysical contribution to ERFACI,sum is larger in the southern NA region than in the northern NA region (71.5 % vs. 54.4 %).

Figures 15 and 16 show the spatial patterns of the contributions and reveal that the contribution from Nd is restricted to the northern part of the domain (including the region of the northern NA subregion) and is generally small elsewhere, except for a small contribution down the west coast of Africa. Somewhat surprisingly, this contribution is not co-located with where the largest changes in Nd are simulated in Fig. 13 but occurs to the north and east. This likely reflects the fact that the cloud fraction is large in this region due to the prevalence of stratocumulus, which results in a large radiative impact from even modest Nd increases. A similar argument can be applied for LWPic changes. The forcing contribution from LWPic changes follows a similar spatial pattern to that from Nd changes and with similar magnitudes, although for LWPic the largest LWPic changes are co-located with the forcing contributions. The Nd contribution is generally dominant in terms of Eq. (4) north of 36 N but also provides a relatively large contribution down the west coast of Africa. The contribution from fc changes is largest in the southern NA region and around the coast of west Africa.

Figure 13Mean percentage increase in Nd, LWPic and fc between PI and PD runs.

Figure 14Time-mean surface ACI forcing (W m−2) calculated using different techniques. (a) Estimated surface ACI forcing calculated using the changes in LWPin-cloud, Nd and cloud fraction between PD and PI simultaneously, according to Eq. (F13); (b) the sum of the estimated contributions from the individual changes (Eq. F15).

3.2.2 Determining the most important meteorological situations for surface ACI forcing

We now attempt to determine the meteorological situations that are most important for aerosol forcing, since this information can then be used to target particular situations for modelling at high resolution or for more detailed comparisons to observations.

Initially, we quantify the forcing as a function of “cloud state”, which we define in terms of the eight possible combinations of low-, mid- and high-altitude clouds (see Fig. 17), which are defined to be consistent with the International Satellite Cloud Climatology Project (ISCCP) definitions (Rossow et al.1999): cloud-top pressure (CTP) > 680 hPa for low-altitude clouds, 680  CTP < 440 hPa for mid-altitude clouds and CTP  440 hPa for high-altitude clouds. For this analysis, clouds are considered present if the cloud fraction is greater than 0.01. The results are presented in the form of 2-D matrix plots showing the contribution to the overall ERFACI for different combinations of PI and PD cloud states (Fig. 18) for both the northern NA and southern NA regions. The contributions take into account the frequency of occurrence of each pairing of PI and PD cloud states so that the values associated with the colours in the plots add up to the overall regional forcing.

Generally, for both regions, the largest negative contributions to ERFACI are associated with low-altitude cloud states (cloud states that have even numbers). To some extent, this is expected since the convection scheme, which is likely one of the sources of higher-altitude clouds, is not currently coded to respond to aerosol and because the lower-altitude clouds are closer to the surface aerosol sources. For the northern NA region, the largest terms are those with the same low cloud state in the PI and PD (i.e. the diagonal terms for the even numbered states in the matrix), although it should be noted that the diagonal terms could involve fc or LWPic changes within these states. The largest overall term occurs when both the PI and PD are in cloud state 8 (i.e. no change in cloud state). This is when low-, mid- and high-altitude clouds are present at the same time. However, if the mid- and high-altitude clouds are thin, it is still possible that the low-altitude cloud has the largest impact here. The next two (joint) largest forcing contributions occur when both the PI and PD have only low clouds (state no. 2) and when low clouds and high clouds are present together (state no. 6).

Figure 15Estimated contributions to the surface ACI forcing from changes in Nd (a), LWPic and fc (b). Figure 14b shows the sum of these three terms.

Figure 16The percentage of the sum of the absolute contributions to the ACI forcing for the terms in Fig. 15. See Eq. (4). Note that the absolute values of the percentages for each grid box add up to 100 %.

However, contributions from transitions between cloud states are not negligible for the northern region. For example, PI to PD transitions between states 1 and 2 (clear sky in the PI and low clouds in the PD) contribute 0.35 W m−2. However, there is a reciprocal positive response of 0.30 W m−2 for transitions between states 2 and 1, which represents the removal of a low cloud that was present in the PI to create clear skies in the PD. There is a certain degree of randomness that is likely in the cloud state changes due to all of the changes not necessarily being driven by aerosol; there is still some model freedom despite the nudging to meteorology since the model is only nudged above the boundary layer and only the winds are nudged. Thus, it is perhaps more useful to consider the net forcing contribution from both the PI to PD transitions and those from the reciprocal transitions, as shown in the bottom row of Fig. 18. Appendix G explains the details of how these are produced. Other transitions have a larger net negative effect, such as the transitions between states 6 and 8 (net effect of 0.18 W m−2). This indicates transitions from low- and high-altitude clouds to low-, mid-, and high-altitude clouds, suggesting that the creation of mid-altitude clouds also has an impact. Nevertheless, these transition terms are considerably smaller than the diagonal (no transition) terms described earlier for the northern NA region.

For the southern NA region, the low-altitude-only state (no. 2) in PI and PD again has a large impact (0.53 W m−2). There is also now a large negative contribution due to transitions from the clear state (no. 1) in the PI to the low-altitude only state 2 in the PD, but again there is a positive forcing due to the reverse situation, with a net contribution of 0.29 W m−2. This transition indicates the creation of additional low-altitude clouds in the PD due to aerosol effects, which is consistent with the finding in the previous section that macrophysical changes to clouds (in particular changes in fc) are more important in this region compared to the northern NA region and suggest that low-altitude clouds are one of the most important cloud types for this. However, as in the northern NA region, pairings with the same cloud states in the PI and PD (the diagonal terms) that involve mid- and high-altitude clouds are also important for forcing, particularly for states 6 and 8, although these are relatively less important than in the northern NA region, suggesting that a higher-altitude cloud has less effect on the forcing in the southern NA region.

3.2.3 Determining the most important cloud types for surface ACI forcing

We now break down the forcing as a function of the PI and PD cloud fractions in order to determine the types of low clouds involved, i.e. whether they are stratocumulus (high cloud fraction) or trade cumulus (low cloud fraction). Given that the contributions involving only cloud state 2 (low-altitude only clouds) were the single largest contribution in Fig. 18 for the southern NA region and provided the joint second-largest contribution for the northern NA (with low clouds also involved for the highest contributor), we examine the forcing contribution for this cloud state only, along with the clear-sky state (no. 1).

For the northern NA region, Fig. 19 shows that by far the largest single term comes from the situation when both the PI and PD are fully overcast (2.24 W m−2). Consistent with the previous figures, this represents mainly the brightening of stratocumulus clouds due to increases in droplet concentrations combined with a smaller macrophysical effect from LWPic increases. There are also some large negative contributions associated with increases in cloud fraction between PI and PD for this region too. The largest net contributions (ranging from 0.1 to 0.37 W m−2) associated with cloud fraction changes come from the creation of fully overcast or fc=0.8 stratocumulus from lower cloud fractions (ranging from 0 to 0.8), indicating that the creation of overcast clouds from more broken stratocumulus and cumulus is playing some role in this region. These contributions sum to 0.91 W m−2.

Figure 17The cloud states used in the following figures. The shading indicates the presence of low-, mid- or high-altitude clouds (see text for the definitions of this) as determined by requiring the model cloud fraction to be larger than 0.01. Cloud state no. 1 is clear sky.


Figure 18Contributions (colours; the highest 10 values in absolute magnitude are also labelled with text) to the surface ACI forcing for the different combinations of cloud state in the PI emission run (x axis) and in the PD emission run (y axis). The cloud states are the eight different possible combinations of low-, mid- and high-altitude clouds (see Fig. 17). The overall contribution is shown, which includes the frequencies of occurrence such that the sum of all of the values gives the overall regional mean contribution. (a) The northern NA region; (b) the southern NA region. The top row shows all of the combinations. In the bottom row plots, the contributions due to PI to PD transitions between states are added to the same transition between PD and PI in order to get the net contribution; see Appendix G for details of this.


For the southern NA region, the contribution from the overcast cloud state in both the PI and PD is much smaller (0.24 W m−2) compared to the northern NA, indicating that pure microphysical (Twomey) changes are less important. The net contribution from zero and overcast state combinations is (-1.50+1.14=-0.36 W m−2) and there are also large net contributions from PI fc values between 0 and 0.4 transitioning to PD values of 0.8–1.0. The creation of fc=0.2 clouds from clear skies is the most important single term associated with cloud fraction changes for the southern NA region, in contrast to the northern NA region. The Twomey effect for fc≤0.2 clouds is also important. This indicates a relatively more important role for trade cumulus or open-cell stratocumulus for the southern NA region. As also suggested from the previous results, this indicates that cloud fraction changes are relatively more important than the Twomey effect for this region compared to the northern NA region.

Figure 19As for Fig. 18 except plotted as a function of the different combinations of areal low cloud fraction of the PI and PD runs.


4 Conclusions

Previous work (Booth et al.2012) has suggested that cloud–aerosol forcing is the dominant driver of multi-decadal SST variability in the North Atlantic region. In this paper, we examined the cloud–aerosol response and surface aerosol forcing of a more modern version of the climate model used in Booth et al. (2012), namely the UKESM1-A, which is the atmosphere-only (i.e. non-ocean-coupled) version of the model used in the latest CMIP intercomparison (CMIP6). We focused on the North Atlantic region and used 1 year of meteorologically nudged model output. The aims were to (i) examine the representation of clouds and cloud–aerosol interactions in order to identify potential biases in cloud properties compared to observations; (ii) quantify the effect of cloud property biases on the SW TOA fluxes; (iii) determine the surface (downward) aerosol forcing response of the model; and (iv) identify the most important meteorological situations/cloud types for the surface forcing. The latter two aims should allow a more targeted evaluation of the model forcing in future work using observations and high-resolution modelling in order to test the magnitude of the forcing from the low-resolution climate model.

4.1 Model evaluation conclusions

The spatial patterns of low-, mid- and high-altitude clouds in the model compare well against observations. However, in the southern North Atlantic (southern NA), the model underestimated the low- and mid-level cloud fractions by 23.9 % and 25.8 %, respectively. Low-altitude cloud biases in the northern NA subregion were positive but lower in magnitude (5.1 %).

Low-altitude cloud fraction biases have the potential to significantly impact aerosol forcing since they are closest in altitude to the aerosol sources. If we assume that the PI cloud cover is biased by a similar amount to the PD cloud cover and assume no cloud adjustments to aerosol, then the bias in forcing will be similar to the bias in fc. The results from this paper suggest that the second assumption is reasonable for the northern NA region because the Twomey effect dominates. Thus, we might expect the 5.1 % low-altitude fc bias there to make a small contribution to any error in forcing bias. Of course, it is also possible that the fc increase between PI and PD is underestimated by the model for this region, but the PI low-altitude fc is too high in order to give an overall PD bias. As such, it is difficult to make firm conclusions.

In the southern NA, fc changes in response to aerosols (adjustments) were large, so the above assumptions are less valid and the effects on forcing even less clear. The negative present-day fc biases could indicate a cloud fraction response to aerosol that is too low, which would cause subsequent negative forcing biases. On the other hand, the too-low fc may mean that the model PI era is in a broken, precipitating cloud regime too often. Such regimes are thought to be more sensitive to aerosols and more prone to produce cloud adjustments (Ackerman et al.2004). In this case, the model forcing values would be too large. The effect of this may be significant given the large bias here (23.9 %). Likewise, the too-large fc values in the northern NA (if also occurring in the PI) might prevent some instances of fc increase between PI and PD and lead to a forcing that is too small, although the overall bias for the northern NA region was only 5.1 %. The presence of clouds can also mask ARI forcing, and hence the fc and LWPic biases might therefore affect the predicted ERFARI magnitude. The larger low- and mid-altitude biases (which are likely to be thicker and hence have a stronger masking effect than high-altitude clouds) in the southern NA combined with the larger ERFARI suggest that this effect would likely be more pronounced there. Here, the fc biases are negative, which would produce a positive ERFARI bias by this mechanism. Mid- and high-altitude clouds can also mask ERFACI forcing from low-altitude clouds (which is likely to be the biggest contributor to forcing). For both the northern NA and southern NA, mid-altitude clouds tend to have negative biases that are larger in magnitude than the positive biases of the high-altitude clouds suggesting the potential for an overall negative ERFACI forcing from this mechanism. However, further work would be needed to quantify the effect of the model fc biases on the aerosol forcing.

In-cloud liquid water path (LWPic) was evaluated vs. the AMSR-E microwave radiometer satellite instrument. The model percentage LWPic biases are small in the northern part of the Atlantic where higher cloud fraction stratocumulus dominates and where the satellite retrievals are likely to be more reliable due to lower rain amounts. In the NE part of the domain (again a region with likely reliable retrievals), there were some positive biases, particularly north of Scandinavia, which indicates that the stratocumulus in this region are too thick in the model. Elsewhere, off the east coast of Florida, there tended to be negative LWPic biases. When the raining data points were filtered out, the biases tended to get worse, but the spatial pattern was preserved, suggesting that the biases are likely real rather than a result of artefacts in the observations or comparison method.

The overestimate of LWPic in the model in the regions dominated by stratocumulus (northern part of the domain) has implications for its ability to simulate the correct cloud effective radius (re) given the correct cloud droplet concentration. This may confuse evaluation efforts using re. An incorrect simulation of cloud droplet size also has implications for the conversion of cloud water into rain (autoconversion), suggesting the model might convert clouds into rain too readily at a given cloud droplet concentration; this will affect rain rates but may also make the response of cloud fraction more sensitive to aerosols and hence enhance aerosol forcing. LWPic is also likely to affect the sensitivity of the cloud albedo to cloud droplet concentration (and hence aerosol), which will also affect the magnitude of aerosol forcing.

Model cloud droplet number concentrations (Nd) were compared to satellite observations from MODIS. The model captures the observed spatial pattern well, suggesting that the processes that govern the removal or dilution of aerosol as it travels eastwards from the North American continent are broadly captured by the model. Such processes are likely to include the scavenging of CCN during precipitation (Wood et al.2017); dilution effects as distance from the sources increases; variations in boundary layer height, which may affect aerosol scavenging due to cloud type and precipitation changes, and may cause the dilution of aerosol concentration over deeper boundary layers; as well as other unidentified meteorological effects. However, the model underestimates Nd over most of the northern NA but overestimates it over the southern NA, which could indicate that there are some issues with aerosol sources and/or scavenging. Cloud processes could also be involved, such as updraft speed distribution errors, problems with the droplet activation scheme or issues relating to droplet removal (evaporation, coagulation, etc.). A caveat to the conclusion that the model has biases in Nd is that there is uncertainty in the MODIS Nd observations that may be of the same magnitude as, or greater than, the model biases (Grosvenor et al.2018b).

The positive model Nd biases associated with the aerosol outflow regions of the US and Europe have the potential to cause a positive bias in the aerosol forcing. In the northern parts of the Atlantic, which are less affected by anthropogenic aerosol, the negative Nd biases may indicate a lack of aerosols from natural sources (or too much scavenging). This would create pre-industrial background aerosol concentrations that are too low, leading to a positive forcing bias (e.g. Carslaw et al.2013). Further investigation into the cause of the spatial Nd pattern using model experiments, an examination of the realism of anthropogenic and natural aerosol sources, and quantification of the MODIS Nd uncertainties (e.g. through comparison with aircraft data for this region) is therefore warranted.

Shortwave top-of-the-atmosphere (SWTOA) radiative fluxes from the model were compared to those from CERES. Again, the model reproduced the observed spatial pattern well, indicating a general good fidelity of overall cloud positions and properties. However, there was also an overestimate in most regions, showing that the clouds were either too bright, occurred too frequently or both.

The main cloud properties that determine SWTOA fluxes are fc, LWPic and Nd. Offline radiative analysis suggested that the fc bias is likely to contribute the most to the SWTOA biases, particularly in the stratocumulus region in the northern Atlantic, suggesting that biases in fc are the most important to address. This result also shows that biases in the response of cloud fraction to aerosol are likely to cause large forcing biases. LWPic biases were determined to be most important in causing the positive SWTOA bias north of Scandinavia and so improvements there are also needed. Nd biases were deemed important in causing SWTOA biases in the northern NA but with only small regions of high impact in the southern NA regions. However, it should be considered that a small impact from biases in a variable on the mean SWTOA flux bias does not preclude a large impact on forcing since the magnitude of the aerosol forcing is a lot smaller (of the order of 1 %–2 %) than the mean SWTOA flux and so small biases can still have the potential to produce a significant forcing error.

4.2 Aerosol forcing conclusions

The ARI and ACI surface aerosol forcings were calculated using the instantaneous model output. In agreement with B12, the magnitude of the ARI forcing over ocean grid boxes in the NA region was generally lower than the ACI forcing. The ACI forcing for the NA region (including land points) is negative with a mean value over the region of 1.7 W m−2, showing that aerosol forcing is likely to be important in terms of the energy received at the ocean surface in this region. The spatial pattern agrees well with that from B12; however, the magnitude is somewhat lower. This likely reflects the steps taken to reduce the aerosol forcing in the UKESM model as described in Mulcahy et al. (2018).

The ACI forcing was decomposed into contributions from changes between PI and PD in Nd, LWPic and fc (ΔNd, ΔLWPic and Δfc, respectively) using an offline calculation of the net surface SW downwelling radiative fluxes. The results showed that ΔNd and ΔLWPic contributed most in the northern part of the NA where overcast stratocumulus clouds are prevalent. Δfc had the highest contribution in the southern subtropical regions where cloud fractions are lower, indicating more broken clouds. We speculate that the higher cloud cover in the northern region allows the ΔNd and ΔLWPic effect to be larger since the associated albedo perturbations can operate over a wider area. In the southern region, it is likely that the broken cloud scenes are more susceptible to aerosol-induced cloud fraction increases. In contrast, in the northern NA region, the overcast clouds cannot increase in cloud fraction much further.

Changes in LWPic and fc can be considered to be cloud macrophysical responses because they are linked to the thermodynamics of liquid water production, whereas changes in Nd can be considered more of a microphysical response since they are mainly caused by aerosol changes, although there is a link between the two. In the northern NA region, the macrophysical changes contribute to approximately 54.4 % of the total ACI forcing, whereas in the southern North Atlantic region they contribute 71.5 %.

Our results are important in the context of Malavelle et al. (2017), who showed that, for an earlier version of this model, the macrophysical responses to volcanic sulfate aerosol perturbations are likely to be small in the impact region of the volcano (north of  50 N). Note that the Malavelle et al. (2017) study quoted responses of the all-sky LWP, which is the product of LWPic and fc. Our Figures 15 and 16 show that with an updated model, using a year of data and considering PI to PD aerosol changes, macrophysical changes (dominated by LWPic changes) have a similar radiative impact to Nd changes for this region. It would be interesting to discover whether the implied larger LWPic response to aerosol in the newer model in the Holuhraun region is now larger than suggested by the observational constraint used in Malavelle et al. (2017). Furthermore, our results show that macrophysical responses are even more important in other NA regions, especially the southern NA. This result highlights that cloud responses to aerosol perturbations in specific regions are not necessarily representative of those elsewhere.

In order to understand the cause of the aerosol-induced LWPic and fc increases in the model (i.e. the macrophysical responses), simulations have been performed where Nd has been prevented from affecting rain formation from cloud liquid through the autoconversion parameterization. Instead, a constant Nd value is used. In these simulations, the changes in LWPic and fc between PI and PD in the NA region drop to nearly zero or even decrease between PI and PD, strongly indicating that the mechanism for the increases in macrophysical quantities in the standard simulations is the suppression of rain. This is consistent with previous studies which have shown that precipitation is a major factor in causing cloud breakup (Stevens et al.1998; Berner et al.2013) mainly due to the stabilization effect of rain evaporation in the lower boundary layer but also with some positive feedback due to the removal of aerosol by rain and the subsequent enhancement of rain due to the lower aerosol concentrations (larger droplets). In a global model, the representation of the thermodynamics of this process will be reliant on the boundary layer scheme, along with input from the cloud microphysics scheme. Since both processes are highly parameterized in global models, it is likely that there is some uncertainty in their representation and thus uncertainty in the response of the clouds to aerosol. Since the contribution to aerosol forcing from changes in macrophysical quantities is very large in this model, this highlights the need for detailed and targeted model evaluation of these processes.

Our results are somewhat different from those obtained in Mülmenstädt et al. (2019) using a similar technique for the ECHAM-HAMMOZ model. For the North Atlantic region, their model produced ERFACI contributions from Nd, LWPic and fc that were located in approximately the same locations in contrast to the fc contribution from our model being predominantly in the southern NA and the Nd and LWPic contributions being further north. This suggests that the lifetime effect (and the associated aerosol-induced precipitation suppression process described above) operates differently between the two models, or that the types and locations of clouds differ between the models. Their model also has very little aerosol forcing in the region north of Scandinavia, whereas our model had quite a large forcing due to Nd and LWPic changes there. It is possible that some of these differences are due to the use of different years or different decomposition techniques. Nevertheless, understanding these model differences may lead to a method of model evaluation to determine which is correct that might help bring down aerosol forcing uncertainty.

4.3 Important meteorological situations: conclusions

Contributions to the surface forcing were calculated for situations composed of each of the eight different combinations (termed here as “cloud states”) of cloudy or clear for low, middle and high cloud altitudes (see Fig. 17). The largest contributions (taking into account the frequency of occurrence as well as the forcing effect) always involved the presence of low-altitude clouds in either the PI or PD suggesting that such clouds are driving most of the aerosol forcing in the model. There was some expectation of this a priori since the convection scheme in the model does not respond to aerosol and it is likely that this would be responsible for a lot of the creation of higher-altitude clouds. However, it is also possible that such clouds are created by the convection scheme and then get handed over to the large-scale cloud scheme and could then be affected by aerosol. It may be that low clouds are more affected simply because they are closer to the surface aerosol sources.

In both regions, the largest forcing contributions came from situations with the same cloud states in the PI and PD. Of these, the largest contributions for the northern NA region (and a considerable contribution for the southern NA region) was when low-, mid- and high-altitude clouds were present, at the same time suggesting that in these regions considerable forcing occurs when higher-altitude clouds are present. Part of the reason for this is likely due to the prevalence of clouds from a range of altitudes in this region. Thus, it may be important to consider the potential radiative effects of mid- and high-altitude clouds such as the shielding of low-altitude clouds, as also indicated in Malavelle et al. (2017). The neglecting of mid- and high-altitude clouds in the radiative calculations presented in this paper may therefore lead to some inaccuracy. On the other hand, the mid- and high-altitude clouds may be sufficiently thin that it has negligible impact. Further work is needed to determine whether this is the case. However, the idea is supported by the fact that the offline calculations of ACI forcing that considered only low-altitude clouds matched the model calculated forcing (that included all cloud types) very closely.

Situations with only low-altitude clouds in the PI and PD were also very important for both the northern NA and southern NA regions, again highlighting the importance of low-altitude clouds for aerosol forcing. Net contributions involving the creation or destruction of a cloud type (low, middle or high clouds) were small for the northern NA region but higher for the southern NA region. For the latter, the creation of low-altitude clouds was implicated in each case. This fits with the result that Δfc was a large contributor to the surface forcing in the southern NA region.

Overall, the results suggest that low-altitude clouds are the largest contributor to aerosol forcing in the model and so improvements to the representation of this and its response to aerosol are likely to yield the biggest improvement in the aerosol forcing estimate in this model. On the other hand, the lack of aerosol awareness of the convection scheme is unrealistic and may lead to missing aerosol–cloud interaction processes in the model. There are some indications that cyclones respond to increased aerosol concentrations by increasing their LWP and hence radiative forcing McCoy et al. (2018). The degree to which the convection scheme is involved in representing the clouds associated with cyclones and hence the likelihood of such aerosol–cloud interactions being missed at climate model resolution could be explored with high-resolution models (e.g. the nested version of the UM) in future work.

4.4 Important cloud types: conclusions

The forcing for grid boxes with low-altitude-only clouds and clear skies was further examined in terms of pairings of cloud fraction from the PI and PD simulations. The aim was to try to understand the PI to PD cloud transitions that occur to produce the model aerosol forcing and therefore gain insight into the types of clouds involved. The results showed that grid boxes with overcast (100 % cloud cover) in both the PI and the PD accounted for by far the largest contribution to forcing in the northern NA region, although there were some smaller contributions from transitions between 60 %–80 % and 100 % cloud cover. This agrees with the result from the forcing contribution calculations that changes in cloud fraction due to aerosol are not particularly important for forcing in this region and suggests that cloud fraction transitions that do occur are mostly likely within the stratocumulus regime. The brightening of stratocumulus clouds by the Twomey effect seems to dominate the forcing in this region.

For the southern NA region, the largest contribution is from transitions from 0 to < 30 % cloud cover, indicating the creation of trade cumulus from clear skies. The creation of 80 %–100 % cloud fraction states from 0 %–20 % states also has an important forcing contribution, suggesting the formation of stratocumulus clouds from clear skies and cumulus to stratocumulus transitions due to the anthropogenic aerosol input. The brightening of overcast clouds is also still very important here but not as important as in the northern NA.

A caveat is that with this analysis is that we only pair the same times in the PI and PD and do not allow for evolution over time. Therefore, the pairings may not be accurate since the Lagrangian trajectories associated with a given pairing may have had different cloud fractions in the time before the snapshot. It may therefore be useful to examine the evolution of cloud fraction, etc. over time using Lagrangian trajectories and examine how this varies with initial aerosol/droplet concentrations (e.g. see Eastman and Wood2016; Eastman et al.2016; Eastman and Wood2018). This could be done in the PD simulation and compared to observations, but also comparisons between the evolution in the PI and PD runs could be made.

4.5 Recommendations

The results from this paper suggest that future studies should target the improvement of different cloud responses to aerosol depending on the location within the North Atlantic. In the northern NA, the Twomey effect was a large contributor to the aerosol forcing. This suggests that constraint of the aerosol forcing using observations will be easier here than in the southern region where cloud adjustments/macrophysical effects were more dominant. Similar approaches to those performed in Malavelle et al. (2017), whereby the model Nd and re responses to known aerosol perturbations (the Holuhrahn volcanic eruption in the case of Malavelle et al.2017) were evaluated using satellite data, may therefore be able to constrain aerosol forcing. Examination of modelled trends in Nd and re over time in response to known aerosol emission changes may also prove useful in this regard. The examination of Nd and re changes alone is advantageous since their natural variability is significantly lower than that of LWPic, fc and SW fluxes (Malavelle et al.2017), making it easier to quantify aerosol-induced signals from the observations. However, our model results also showed a fairly large influence on aerosol forcing from increases in LWPic (i.e. cloud thickening) in the northern NA region, which may complicate such efforts.

In the southern NA, the accuracy of the cloud fraction response in the model should be evaluated, with a particular focus on low clouds and the creation of trade cumulus and stratocumulus. As mentioned earlier, this is likely to be more difficult than evaluating the Twomey effect due to the larger number of processes involved. Ideally, the individual processes would be targeted for model evaluation with the formation of rain via the autoconversion process a prime first target since this was shown to be the cause of the fc and LWPic response to aerosol in our model. Ways forward with this might include the use of observational data from the ground-based ARM site on the Azores and data from aircraft observations that have also taken place there. This is ideally located since it is near a location where the older version of the model predicts a large contribution to the aerosol forcing from cloud fraction changes but the newer version less so. Thus, the observations might help to determine which is the most realistic. Combined with satellite data, the long time series observations of aerosols, Nd, LWPic, fc and rain rates can be used to evaluate the relationships between these variables. Separating a causative signal due to the aerosol alone from that due to meteorology, etc. is difficult, but it may be possible to use the techniques described in McCoy et al. (2020). In that paper, several meteorological drivers are controlled for via binning and multi linear regression, along with using the model to estimate and subtract the influence of non-aerosol-induced changes. However, it is possible that even an evaluation of the observed relationship between rain rates and other variables, without separating the causative effects of aerosol, will prove useful in constraining the cloud adjustments.

A complementary approach is to model this region using a high-resolution nested version of the model, which would also allow for the use of the more sophisticated CASIM (Cloud AeroSol Interaction Microphysics) microphysics scheme (e.g. see Grosvenor et al.2017; Stevens et al.2018; Miltenberger et al.2018a, b; Gordon et al.2018). The assumption would be that the cloud responses to aerosol of this would be more accurate than the global model resolution simulations, although its performance should be tested using the observations.

In order to avoid issues regarding the representativeness of a limited number of years of meteorology, it would be useful to do future work such as that performed here using multiple years of data. Noise in some of the fields shown here may also be reduced by employing an ensemble approach as in Liu et al. (2018).

Finally, work to evaluate and improve the accuracy of the satellite evaluation in this paper would be very useful. For example, to deal with the issues of evaluating the model LWP in situations with precipitation, a microwave radiometer simulator could be implemented into the model, which would be able to estimate the total attenuation of the 37 GHz channel (used by AMSR-E to retrieve LWP), taking into account the differing attenuation strengths of cloud water and rainwater (Lebsock and Su2014). This attenuation could then be directly compared to that from AMSR-E. Comparing the global model evaluation results with and without the convective contribution to those from a high-resolution model (where the convective TLWP will be explicitly resolved) might help to assess the role of the convective parameterization on the evaluation of model TLWP. This could be combined with the microwave radiometer simulator mentioned above to help overcome the issue of retrieval problems in raining conditions too.

Much greater confidence in Nd retrievals would gained through further validation using in situ aircraft data. Cloud fraction could be evaluated with additional cloud fraction satellite data sets and ground-based data to improve the confidence in the result. More sophisticated analyses such as 2-D histograms of cloud fraction and LWPic would assess this important cloud fraction variable at different cloud thicknesses and may help to isolate issues in particular regimes. Model improvements to the subgrid cloud scheme might be considered, such as a link between the boundary layer turbulence and the width of the subgrid relative humidity distribution. Comparisons to high-resolution models may also help with this.

Appendix A: Details on the satellite data sets

A1 Droplet concentration

Nd can be estimated using satellite retrievals of τc and re (Han et al.1998; Brenguier et al.2000; Boers et al.2006; Bennartz2007; Grosvenor et al.2018b) made using observations from the visible and shortwave infrared wavelengths from instruments like MODIS (Nakajima and King1990). Here, we use a data set based on the methods from Grosvenor et al. (2018a), which represented some modifications to the methods described in Grosvenor and Wood (2014). The methodology used here differs slightly from that used in Grosvenor et al. (2018a) in that data are not filtered to remove data points with τc of less than 5 and the correction for the vertical penetration depth bias proposed in Grosvenor et al. (2018a) is not applied. The 3.7µm re is used for the Nd calculations, which has been suggested to be less prone to errors due to cloud heterogeneity (Zhang et al.2012; Z. Zhang et al.2016; Grosvenor et al.2018b). The data set excludes 1×1 data points with mean solar zenith angle (SZA) greater than 65, mean cloud top heights greater than 3.2 km, liquid cloud fractions less than 80 % and for which the maximum sea-ice areal coverage over a moving 2-week window exceeded 0.001 %. The sea-ice data used were the daily 1×1 version of the “Sea Ice Concentrations from Nimbus-7 SMMR and DMSP SSM/I-SSMIS Passive Microwave Data, Version 1” data set (Cavalieri et al.1996). The Nd data set (without sea-ice screening) is available; see Grosvenor and Wood (2018).

Appendix B: Bias metrics

Following Gustafson and Yu (2012), the NMBF is defined as

(B1) 100 × M O - 1 , if M O and  M M = O O 100 × 1 - O M , if M < O and M M = O O Undefined, if M M O O , i.e. the signs of M and O differ ,

where M is the mean of the model values and O is the mean of the observed values. In this paper, this refers to the spatial mean time-averaged values. The NMAEF is defined as

(B2) 100 × M i - O i O i , if M O and M M = O O 100 × M i - O i M i , if M < O and M M = O O Undefined , if M M O O ,

where Mi are the time-averaged model values and Oi the time-averaged observed values for spatial location i.

Appendix C: Aerosol properties

Figure C1 shows an evaluation of the 550 nm model aerosol optical depth (AOD) using data from MODIS. The MODIS data are a combination of the 550 nm Dark Target and Deep Blue product “Dark_Target_Deep_Blue_Optical_Depth_550_Combined” from the level-3 product (Levy et al.2013). The model shows positive biases in equatorial Africa but negative biases in northern Africa and in the ocean to the west of there. This indicates some issues with dust in the model. Further north, there are small positive biases in the aerosol outflow region to the east of the US. Model values are also too high in the region to the north and west of the UK, and over the Mediterranean region. This may indicate an overestimate of aerosol mass from the land sources of pollution or a lack of scavenging. However, there are potential issues with the comparison between MODIS and the model in these regions since MODIS retrievals are only possible in cloud-free regions and model values are taken in both cloudy and cloud-free regions. Since the cloud coverage of this region is very high (see Figs. 1, 2 and 3), this may cause sampling biases, and further work is needed to examine the effects of this. Overall, model NMBF error values are 8.5 %, 5.6 % and 18.3 % over the NA, northern NA and southern NA regions, respectively.

Figure C1Time-mean 550 nm aerosol optical depth (AOD) model evaluation. The MODIS data are from daytime overpasses of the Aqua satellite using the combined Dark Target and Deep Blue level-3 product. Model local times within 3 h on either side of 13:30 LT are used to approximately match the satellite overpass times.

We now examine the changes between PI and PD for various aerosol-related model fields. Figures C2 and 13 show that the spatial patterns of the changes in AOD, CCN at 0.2 % supersaturation using aerosol from the lowest model layer (CCN0.2 %) and Nd are very similar, suggesting that changes in both AOD and CCN0.2 % are a good proxy for Nd changes in this region. In terms of aerosol composition (Fig. C3), the largest changes in column aerosol mass (the vertically integrated aerosol mass concentrations for all size modes) occur for sulfate and BC. Sulfate changes occur further north than BC changes, with the latter occurring mostly over southern Europe and northern Africa. Sulfate changes dominate over the ocean except at around 18 N where there is a band of larger BC change likely associated with outflow from Africa. Similarly to BC, OM increases also occur over Africa and over the Atlantic Ocean at 18 N, with both BC and OM potentially contributing to the increase in CCN0.2 % and Nd over the Atlantic that region. Sulfate changes likely dominate the Nd changes further north in the Atlantic. OM aerosol mass decreases in the northern part of the Atlantic. Dust and sea-salt column mass changes are very small in comparison to the other aerosol types for this region.

Figure C2Mean percentage increase (between PI and PD model runs) in (a) 550 nm AOD and (b) CCN at 0.2 % supersaturation using aerosol from the lowest model layer.

Figure C3As for Fig. C2 except for changes in (a) column sulfate aerosol mass, (b) column BC aerosol mass, (c) column OM aerosol mass, (d) column dust mass and (e) column sea-salt aerosol mass.

Appendix D: Evaluation of grid-box mean LWP

Here, we show plots similar to Figs. 5 and 6 but for the all-sky (cloudy and clear contributions) LWP, i.e. without applying the step of dividing by the cloud fraction. Table D1 summarizes the results for this and the other plots in this section. The bias pattern is very similar to that of the in-cloud values, suggesting that the conversion between LWP and LWPic (by dividing by the low-altitude COSP-CALIPSO cloud fraction) does not greatly affect the evaluation, probably because the low cloud fraction biases are generally very small (Fig. 1). However, for the regions of negative bias between the Equator and 18 N, the biases are larger for LWP than for LWPic because there is a negative low cloud fraction bias here, which acts to increase the model LWPic values relative to the observed ones to give a lower LWPic bias. When filtering to include only data points with fLWP>0.99 (Fig. D2), the bias pattern is again similar to that for LWPic.

As discussed earlier, AMSR-E observations are potentially biased when it is raining. This is because assumptions are made about the partitioning between LWP and RWP in order to facilitate the retrieval since rainwater attenuates the microwave signal more strongly than liquid droplets (Lebsock and Su2014). It is assumed that rainwater does not occur for water paths below 180 g m−2, which may lead to inaccuracies since the true partitioning value is likely to vary. Because of these issues, some previous model evaluation studies (e.g. Furtado et al.2016) have chosen to compare the TLWP (the sum of LWP and RWP) to the TLWP provided in some microwave-based products. We also do this in Fig. D3; for this plot, we only use the large-scale RWP and not the convective RWP (or LWP). However, it should be borne in mind that if the assumed partitioning threshold is incorrect, the TLWP value retrieved by the satellite will also be incorrect, and so it is unclear whether this leads to a more accurate model-to-satellite evaluation.

The results show that the addition of RWP (compare Figs. D1 and D3) enhances both the model and the observations, so that the pattern of bias remains similar. This suggests that it does not matter greatly whether LWP or TLWP are used for model evaluation. The magnitudes of the negative model biases for TLWP are slightly enhanced and those of the positive biases reduced relative to those for LWP, though. The NMBFs for the northern NA region are 16.4 % for TLWP and 8.3 % for LWP. For the southern NA region, the corresponding values are 65.7 % and 39.5 %.

Figure D4 shows the ratio of the TLWP from the model's convection parameterization to that of the large-scale plus convective TLWP. In most of the southern part of the NA region, the convective TLWP accounts for the majority of the TLWP. However, it is not straightforward to decide whether the convective TLWP from the model is physically meaningful and whether it should be included when comparing to the microwave instruments. On the one hand, LWP and RWP from convection in the real world will be detected by the instruments, but on the other hand it is unclear whether the condensed water from the convection parameterization in models properly represents this. Liquid water content and vapour are detrained to the environment from the convection scheme and incorporated into the large-scale cloud scheme (see UM Documentation Paper 030; hereafter UMDP030;, last access: 11 December 2020). It may therefore be the case that the large-scale cloud amount is the more appropriate quantity even if there is condensate associated with the convection scheme and that double counting would occur if also using the convective values (UMDP30). However, there is also some convective condensate from the convective core that is not transferred to the large-scale scheme. Another point of note is that LWP associated with the convection scheme is only used by the radiation scheme for shallow convective clouds (clouds with geometrical depths less than 500 m over land and 1500 m over ocean), meaning that the majority of liquid water in deep clouds has no effect on radiation. Given the uncertainty, here we examine the effect of adding the convective LWP on the model evaluation (Fig. D5). As might be expected from Fig. D4, its inclusion leads to much larger model values, particularly in the southern parts of the region. The bias pattern also changes, so that large positive biases occur almost everywhere, compared to negative biases in the south and positive biases in the north when using only large-scale TLWP (Fig. D3). Given this, it would be useful for future studies to determine whether any of the convective TLWP should be included; ideas for how to make progress on this are discussed in Sect. 4.5.

The effect of the further addition of convective RWP is shown in Fig. D6. This increases the model values somewhat in the convective regions but not by a large relative amount. As might be expected, the model biases are also therefore increased but not to the same extent as the addition of the convective LWP; the spatial pattern of model bias remains very similar to that in Fig. D5.

We also test the effect of filtering using fLWP>0.99 for the model and satellite data when including the large-scale RWP and the convective LWP and RWP (Fig. D7). The filtering reduces the large positive biases in the southern part of the domain (which were mainly due to the addition of the convective LWP) considerably, presumably because most of the grid points with convection are removed. The spatial patterns of the biases here are similar to those in Fig. D2, where only the LWP from the large-scale cloud scheme was used, but the biases in the southern part of the domain have changed from being slightly negative to slightly positive, while the biases in the northern part have become a little more positive. These relatively small changes suggest that whether the convective LWP, which was the biggest contributor to the total water path, is used is not as critical when filtering for non-precipitating situations. However, there are still large overall differences in the NMBF values between Figs. D2 and D7 (see Table D1) such that whether RWP or convective LWP/RWP is used or not still has significant effects on model evaluation attempts and therefore should be an avenue of future investigation.

Figure D1Time-mean all-sky (i.e. including cloudy and clear regions) LWP model evaluation for both day and night overpasses. Panel (c) includes a contour of the 0.9 value of fLWP; see Fig. 4 for the full map of fLWP for reference.

Figure D2As for Fig. D1 except both the model and satellite data have been filtered before time averaging to only include data points for which fLWP is greater than 0.99. This quantity is denoted as LWP0.99.

Figure D3As for Fig. D1 except with the addition of the all-sky rainwater path (RWP) for both the model and AMSR-E satellite data. No filtering for fLWP has been applied.

Figure D4Ratio of the TLWP (where TLWP is the sum of LWP and RWP) from the model convection scheme to the total (from the convection and large-scale scheme) TLWP. The RWP from the convection scheme is calculated from the convective rain rate diagnostic by assuming a raindrop size distribution and fall speed relationship (see Furtado et al.2016, for details).

Figure D5As for Fig. D3 except with the addition of the all-sky LWP from the convection parameterization for the model.

Figure D6As for Fig. D5 except with the addition of the all-sky RWP from the convection parameterization of the model.

Figure D7As for Fig. D6 except that the model and AMSR-E data have been filtered before time averaging to only include data points for which fLWP is greater than 0.99 (as in Fig. D2).

Table D1Model evaluation statistics for the various subregions using time-averaged data. See Gustafson and Yu (2012) for details of the NMBF and the NMAEF. r is the spatial correlation coefficient between the model and observed time averages. All values are area weighted to account for the variation in area of the model grid boxes.

Download Print Version | Download XLSX

Appendix E: Removing aerosol impact on rain autoconversion

Figures E1 and E2 show the percentage increases between PI and PD for LWPic and fc (ΔLWPic and Δfc), respectively, when aerosols are prevented from affecting the rain autoconversion process. This is done by setting Nd in the autoconversion process equation to a constant value of 300 cm−3 over land and 100 cm−3 over oceans in both the PI and PD runs. These are termed the ConstantNdAutoCon runs here. Most oceanic regions showed positive ΔLWPic and Δfc in the full model and near-zero or negative changes for ConstantNdAutoCon. Table 3 shows that between the standard and ConstantNdAutoCon runs ΔLWPic reduced from 3.2 % to 5.0 % for the NA region and from 2.7 % to 0.43 % for the northern NA. For the southern NA region, there was a small negative ΔLWPic in the standard runs (0.59 %) and a very similar value for ConstantNdAutoCon (0.50 %) consistent with the idea that aerosols have little impact on LWPic in this region, as discussed earlier. Respective Δfc values for the standard and ConstantNdAutoCon runs were 1.5 % and 0.14 % for the NA region; 0.94 % and 0.02 % for the northern NA; and 2.7 % and 0.75 % for the southern NA. These results suggest that most of the PI to PD increases in the macrophysical cloud properties (LWPic and fc) were due to the impact of aerosols on the rain autoconversion process, likely via the precipitation suppression effect of enhanced aerosol.

However, in the region of the Atlantic near the Equator, Fig. E2 suggests that the increase in fc between the PI and PD runs was not due to the impact of aerosol on the autoconversion process since there are still positive Δfc values for the ConstantNdAutoCon runs. We hypothesize that, as far as allowed by the nudging (wind-only nudging applied above the boundary layer), the fc increases in this region may be related to other aerosol impacts such as the ARI, ACI or semi-direct effects, via thermodynamical or dynamical changes.

Figure E1Mean percentage increase in LWPic between PI and PD runs for (a) the full run and (b) the run where aerosol has been prevented from affecting the rain autoconversion.

Figure E2As for Fig. E1 except for fc changes.

Appendix F: Offline shortwave flux calculation and partitioning

Here, we give details of the offline calculations used to estimate the SW fluxes, which are needed to estimate the contributions to the forcing from the individual changes in cloud properties between the PI and PD simulations.

An estimate of the cloud optical depth (τc) can be made following Eqs. (1) and (5) of Grosvenor et al. (2018b):

(F1) τ c = z base z top 3 Q ext 4 ρ w L ( z ) r e ( z ) d z ,

where Qext is the scattering efficiency, which we assume to have a constant value of 2; this has been shown to be the case for droplet radii that are much larger than the wavelength of light concerned (Bennartz2007). ρw is the density of liquid water, L(z) is the liquid water content, re(z) is the effective radius, z is height, zbase is cloud base height, and ztop is cloud top height.

We assume that the clouds are adiabatic (or some constant fraction of adiabatic) so that their liquid water increases linearly with height, and it is assumed that Nd is constant throughout their depth. Observations suggest that both are valid assumptions for stratocumulus clouds (Albrecht et al.1990; Zuidema et al.2005; Painemal and Zuidema2011; Miles et al.2000; Wood2005). The adiabatic assumption means that

(F2) L ( z ) = f ad c w z ,

where fad is the adiabatic fraction, which is assumed constant with height and with a value of 1.0. cw(T,P) is the rate of increase of liquid water content with height (dL∕dz, with units kg m−4) and is referred to as the “condensation rate” in Bennartz (2007) or the “water content lapse rate” in Painemal and Zuidema (2011). See Ahmad et al. (2013) for a definition. A constant temperature of 278 K is used for the temperature (T) in the calculation of cw, along with a constant pressure (P) of 850 hPa. This is an approximation since these values would vary depending upon cloud height and location. However, since we only consider low clouds here and because in Eq. F1 the dependence of τc upon cw is very weak (τ1/cw1/6), this introduces very little error.

In the radiative scheme of the UKESM model the parameterization of Liu et al. (2008) is used, which makes the width of the droplet size distribution (assumed to be represented by a gamma function) a function of Nd and L. For consistency with the model, we also apply this parameterization to our offline radiative calculations. From Mulcahy et al. (2018),

(F3) r e ( z ) = β m ( z ) 3 L ( z ) 4 π ρ w N d 1 3 ,

where βm is a parameter related to the droplet distribution width, which is parameterized in Liu et al. (2008) as

(F4) β m ( z ) = x N d L ( z ) y .

x and y are constants set at 0.0266 (with units kgy, which have been converted from the value of 0.07 gy from Liu et al. (2008)) and 0.14, respectively, following Liu et al. (2008) and Mulcahy et al. (2018). NB: βm is related to the k parameter, which is often used to describe the droplet distribution width (e.g. Martin et al.1994), as βm=1/k(1/3).

Using Eq. (F3) to substitute for re in Eq. (F1) and including Eqs. (F2) and (F4) gives

(F5) τ c = B z base z top z y + 2 3 d z ,


(F6) B = 3 Q ext f ad c w y + 2 3 4 ρ w x 3 4 π ρ w 1 / 3 N d y - 1 3 .

Integrating over height gives

(F7) τ c = B H y + 5 3 y + 5 3 ,

where H is the depth of the cloud (ztopzbase). H can be determined from the in-cloud LWP under the adiabatic assumption:

(F8) H = 2 LWP ic f ad c w 1 2 ,

which follows from integrating Eq. (F2) over the depth of the cloud. This now allows τc to be calculated from the cloud LWPic and Nd.

The cloud albedo (Ac) is then estimated using Eq. (24.38) of Seinfeld and Pandis (2006), which is based on the two-stream approximation for a non-absorbing, horizontally homogeneous cloud:

(F9) A c = τ c τ c + 7.7 .

The shortwave downwards flux at the surface (SWdownSURF) for a given cloud fraction (fc) can then be estimated from the cloudy and clear-sky fluxes (SWdownSURFcloudy and SWdownSURFclear, respectively) using

(F10) SW downSURF = f c SW downSURFcloudy + 1 - f c SW downSURFclear .

SWdownSURFclear is estimated from the incoming TOA SW flux (SWdownTOA) using

(F11) SW downSURFclear = T atmos SW downTOA ,

where we have assumed a constant clear-sky transmissivity (Tatmos). The cloudy-sky surface flux is calculated by assuming that the flux reaching the cloud top is equal to SWdownSURFclear:

(F12) SW downSURFcloudy = T atmos SW downTOA 1 - A c .

Thus, we now have a function SWdownSURF (Nd, fc, LWPic) that estimates the surface downwelling SW flux as a function of the three cloud variables of interest. This allows us to estimate the surface ERFACI following Eqs. (2) and (1) as

(F13) ERF ACI , all = ERF PD - ERF PI = SW PD clean + cloudy - SW PD clean + clear - SW PI clean + cloudy + SW PI clean + clear = SW PD clean + cloudy - SW PI clean + cloudy = SW downSURF ( N d PD , f c PD , LWP ic PD ) - SW downSURF ( N d PI , f c PI , LWP ic PI ) ,

where we have assumed that the PI and PD SWclean+clear values are the same.

The forcing contributions from the changes to the individual cloud parameters are then estimated using

(F14) ERF N d = 0.5 ( ERF N d , PI base + ERF N d , PD base ) ERF f c = 0.5 ( ERF f c , PI base + ERF f c , PD base ) ERF LWP ic = 0.5 ( ERF LWP ic , PI base + ERF LWP ic , PD base ) ,


(F15) ERF N d , PI base = SW downSURF ( N d PD , f c PI , LWP ic PI ) - SW downSURF ( N dPI , f c PI , LWP ic PI ) ERF f c , PI base = SW downSURF ( N d PI , f c PD , LWP ic PI ) - SW downSURF ( N d PI , f c PI , LWP ic PI ) ERF LWP ic , PI base = SW downSURF ( N d PI , f c PI , LWP ic PD ) - SW downSURF ( N d PI , f c PI , LWP ic PI ) .

Here, all of the PI cloud property values have been used as a baseline value but with the PI value for either Nd, fc or LWPic replaced with the PD value. A similar calculation is done using the PD values as baselines and replacing with one of the PI values:

(F16) ERF N d , PD base = SW downSURF ( N d PD , f c PD , LWP ic PD ) - SW downSURF ( N d PI , f c PD , LWP ic PD ) ERF f c , PD base = SW downSURF ( N d PD , f c PD , LWP ic PD ) - SW downSURF ( N d PD , f c PI , LWP ic PD ) ERF LWP ic , PD base = SW downSURF ( N d PD , f c PD , LWP ic PD ) - SW downSURF ( N d PD , f c PD , LWP ic PI ) .

This follows the work of Mülmenstädt et al. (2019), who found that an average of values obtained from using both the PI and PD values as baselines was more accurate than only using say the PI as a baseline. If fc is zero in the baseline state (i.e. PI for Eq. F15, PD for Eq. F16), then Nd and LWPic are undefined. Therefore, in order to calculate the effect of an increased fc between the baseline and perturbed state (i.e. PD for Eq. F15, PI for Eq. F16), the Nd and LWPic values from the perturbed state are used in such cases.

Appendix G: Net contributions to forcing from cloud state and cloud fraction state transitions

Section 3.2.2 and 3.2.3 showed the contributions to the surface forcing from different combinations of PI and PD cloud states (low-, mid- and high-altitude cloud combinations) and cloud fraction. We explained that it is useful to consider the net forcing contribution from both the PI to PD transitions and those from the reciprocal transitions. That is, if gi and gj are two different cloud states (where the possible states run from 1 to 8; see Fig. 17), the net contribution (ERFcontNET) for that pairing of cloud states is given by

(G1) ERFcont NET g i , g j = ERFcont g i PI , g j PD + ERFcont g j PI , g i PD g j > g i ,

where ERFcontgiPI,gjPD is the forcing contribution from pairings of gi in the PI simulation and gj in the PD simulation, and ERFcontgjPI,giPD is the contribution from pairings of gj in the PI simulation and gi in the PD simulation.

A similar equation can be formulated using the five different cloud fraction bins instead of cloud states.

Appendix H: Results from an alternate year

The results in the main body of the paper were based on meteorology and emissions from the period of 28 March 2009 to 28 March 2010. It is possible that the results presented vary depending on the chosen year since meteorology, cloud fields, etc. vary from year to year. To address this, we have run the PI and PD simulations for an additional year and present results here for the period of 28 March 2010 to 28 March 2011. For brevity, we only present selected results focusing on those that are most likely to be subject to variability and noise issues, which are (1) the ARI and ACI aerosol forcings (Fig. 12); (2) the decomposition of the aerosol forcing into contributions from changes in Nd, LWPic and fc between the PI and PD (e.g. Fig. 15); (3) the contributions from changes between different cloud states and between cloud fractions (Figs. 18 and 19).

Figure H1 shows that the pattern and magnitude of both ERFARI and ERFACI are very similar for the alternate year with only slight differences (compare to Fig. 12): the region of ERFARI forcing off the coast of the US extends further east across the Atlantic; the main region of ERFACI around Newfoundland is still present but extends further south and less to the southeast; the region of ERFACI north of Scandinavia has a similar spatial pattern but is enhanced somewhat in the alternate year.

Figure H1As for Fig. 12 except for the 2010–2011 period.

Figure H2 (compare to Fig. 15) shows the contributions to ERFACI from the changes in Nd, LWPic and fc between the PI and PD for the alternate year. It shows that there are some small differences in the spatial patterns, but the overall patterns of the contributions remain very similar with Nd and LWPic changes dominating over the northern NA region and fc changes dominating over the southern NA as was the case for the original chosen time period.

A figure equivalent to Fig. 18 for the alternative year showing the contributions from changes between different cloud states (not shown) reveals a very similar pattern with very similar magnitudes of contribution. Likewise, the alternative version of Fig. 19 (contributions from changes between cloud fractions; not shown) is again very similar to original version in terms of both the pattern and magnitude of the contributions.

Figure H2As for Fig. 15 except for the 2010–2011 period.

Data availability

Raw model data are kept on tape archive available through JASMIN (, last access: 11 December 2020) service. Please see, last access: 11 December 2020 for details on how to arrange access to Met Office data via JASMIN.

Author contributions

DPG designed and ran the simulations, analysed the output, processed the satellite data and produced the text and figures. KSC helped to analyse the model output and provided feedback and edits to manuscript drafts.

Competing interests

The authors declare that they have no conflict of interest.

Special issue statement

This article is part of the special issue “Marine aerosols, trace gases, and clouds over the North Atlantic (ACP/AMT inter-journal SI)”. It is not associated with a conference.


We acknowledge use of the MONSooN system, a collaborative facility supplied under the Joint Weather and Climate Research Programme, a strategic partnership between the Met Office and the Natural Environment Research Council. Version 3 of the CALIPSO-GOCCP (GCM Oriented Cloud Calipso Product; Chepfer et al.2010) was obtained from AMSR-E LWP data are produced by Remote Sensing Systems and sponsored by the NASA Earth Science MEaSUREs DISCOVER Project and the NASA AMSR-E Science Team. Data are available at, last access: 11 December 2020. The MODIS data were obtained from NASA's Level 1 and Atmosphere Archive and Distribution System (LAADS) (, last access: 11 December 2020). CERES-EBAF data were taken from and are the monthly averaged product of observed TOA for which the TOA net flux is constrained to the ocean heat storage.

Finally, we thank the two anonymous referees who helped to improve the manuscript through their constructive comments.

Financial support

This research has been supported by the National Environment Research Council (North Atlantic Climate System: integrated Study (ACSIS) project). Daniel P. Grosvenor was supported by the National Environmental Research Coun40 cil (NERC) national capability grant for the The North Atlantic Climate System Integrated Study (ACSIS) programme (grant no. NE/N018001/1) via NCAS.

Review statement

This paper was edited by Hailong Wang and reviewed by two anonymous referees.


Abdul-Razzak, H. and Ghan, S. J.: A Parameterization of Aerosol Activation: 2. Multiple Aerosol Types, J. Geophys. Res., 105, 6837–6844,, 2000. a

Abel, S. J. and Shipway, B. J.: A comparison of cloud-resolving model simulations of trade wind cumulus with aircraft observations taken during RICO, Q. J. R. Meteorol. Soc., 133, 781–794,, 2007. a

Ackerley, D., Booth, B. B. B., Knight, S. H. E., Highwood, E. J., Frame, D. J., Allen, M. R., and Rowell, D. P.: Sensitivity of Twentieth-Century Sahel Rainfall to Sulfate Aerosol and CO2Forcing, J. Climate, 24, 4999–5014,, 2011. a

Ackerman, A. S., Kirkpatrick, M. P., Stevens, D. E., and Toon, O. B.: The impact of humidity above stratiform clouds on indirect aerosol climate forcing, Nature, 432, 1014–1017,, 2004. a, b, c

Ahmad, I., Mielonen, T., Grosvenor, D. P., Portin, H. J., Arola, A., Mikkonen, S., Kühn, T., Leskinen, A., Joutsensaari, J., Komppula, M., Lehtinen, K. E. J., Laaksonen, A., and Romakkaniemi, S.: Long-term measurements of cloud droplet concentrations and aerosol-cloud interactions in continental boundary layer clouds, Tellus B, 65, 20138,, 2013. a

Albrecht, B. A.: Aerosols, Cloud Microphysics, and Fractional Cloudiness, Science, 245, 1227–1230,, 1989. a

Albrecht, B. A., Fairall, C. W., Thomson, D. W., White, A. B., Snider, J. B., and Schubert, W. H.: Surface-Based Remote-Sensing Of The Observed And The Adiabatic Liquid Water-Content Of Stratocumulus Clouds, Geophys. Res. Lett., 17, 89–92,, 1990. a

Andreae, M. O., Jones, C. D., and Cox, P. M.: Strong present-day aerosol cooling implies a hot future, Nature, 435, 1187–1190,, 2005. a

Archibald, A. T., O'Connor, F. M., Abraham, N. L., Archer-Nicholls, S., Chipperfield, M. P., Dalvi, M., Folberth, G. A., Dennison, F., Dhomse, S. S., Griffiths, P. T., Hardacre, C., Hewitt, A. J., Hill, R. S., Johnson, C. E., Keeble, J., Köhler, M. O., Morgenstern, O., Mulcahy, J. P., Ordóñez, C., Pope, R. J., Rumbold, S. T., Russo, M. R., Savage, N. H., Sellar, A., Stringer, M., Turnock, S. T., Wild, O., and Zeng, G.: Description and evaluation of the UKCA stratosphere–troposphere chemistry scheme (StratTrop vn 1.0) implemented in UKESM1, Geosci. Model Dev., 13, 1223–1266,, 2020. a

Ba, J., Keenlyside, N. S., Latif, M., Park, W., Ding, H., Lohmann, K., Mignot, J., Menary, M., Wouters, O. H. O. B., Melia, D. S., Oka, A., Bellucci, A., and Volodin, E.: A multi-model comparison of Atlantic multidecadal variability, Clim. Dynam., 43, 2333–2348,, 2014. a

Bennartz, R.: Global assessment of marine boundary layer cloud droplet number concentration from satellite, J. Geophys. Res.-Atmos., 112, 1–16,, 2007. a, b, c

Berner, A. H., Bretherton, C. S., Wood, R., and Muhlbauer, A.: Marine boundary layer cloud regimes and POC formation in a CRM coupled to a bulk aerosol scheme, Atmos. Chem. Phys., 13, 12549–12572,, 2013. a, b

Bodas-Salcedo, A., Webb, M. J., Bony, S., Chepfer, H., Dufresne, J.-L., Klein, S. A., Zhang, Y., Marchand, R., Haynes, J. M., Pincus, R., and John, V. O.: COSP: Satellite simulation software for model assessment, B. Am. Meteorol. Soc., 92, 1023–1043,, 2011. a

Boers, R., Acarreta, J. R., and Gras, J. L.: Satellite monitoring of the first indirect aerosol effect: Retrieval of the droplet concentration of water clouds, J. Geophys. Res., 111, D22208,, 2006. a

Booth, B. B. B., Dunstone, N. J., Halloran, P. R., Andrews, T., and Bellouin, N.: Aerosols implicated as a prime driver of twentieth-century North Atlantic climate variability, Nature, 484, 228–232,, 2012. a, b, c

Boucher, O., Randall, D., Artaxo, P., Bretherton, C., Feingold, G., Forster, P., Kerminen, V.-M., Kondo, Y., Liao, H., Lohmann, U., Rasch, P., Satheesh, S. K., Sherwood, S., Stevens, B., and Zhang, X. Y.: Clouds and Aerosols, in: Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Stocker, T. F., Qin, D., Plattner, G.-K., Tignor, M., Allen, S. K., Boschung, J., Nauels, A., Xia, Y., Bex, V., and Midgley, P. M., Cambridge University Press, Cambridge, UK and New York, USA, 571–658, 2013. a

Boutle, I. A., Abel, S. J., Hill, P. G., and Morcrette, C. J.: Spatial variability of liquid cloud and rain: Observations and microphysical effects, Q. J. R. Meteorol. Soc., 140, 583–594,, 2014. a

Brenguier, J.-L., Pawlowska, H., Schuller, L., Preusker, R., Fischer, J., and Fouquart, Y.: Radiative properties of boundary layer clouds: Droplet effective radius versus number concentration, J. Atmos. Sci., 57, 803–821,<0803:RPOBLC>2.0.CO;2, 2000. a

Bretherton, C. S., Blossey, P. N., and Uchida, J.: Cloud droplet sedimentation, entrainment efficiency, and subtropical stratocumulus albedo, Geophys. Res. Lett., 34, L03813,, 2007. a, b

Brown, A. R., Beare, R. J., Edwards, J. M., Lock, A. P., Keogh, S. J., Milton, S. F., and Walters, D. N.: Upgrades to the boundary-layer scheme in the met office numerical weather prediction model, Boundary-Layer Meteorol., 128, 117–132,, 2008. a

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. a

Cavalieri, D. J., Parkinson, C. L., Gloersen, P., and Zwally, H. J.: Sea Ice Concentrations from Nimbus-7 SMMR and DMSP SSM/I-SSMIS Passive Microwave Data, Version 1, NASA National Snow and Ice Data Center Distributed Active Archive Center, Boulder, Colorado USA,, 1996. a

Chepfer, H., Bony, S., Winker, D., Cesana, G., Dufresne, J. L., Minnis, P., Stubenrauch, C. J., and Zeng, S.: The GCM-oriented CALIPSO cloud product (CALIPSO-GOCCP), J. Geophys. Res.-Atmos., 115, 1–13,, 2010. a, b

Cox, P. M.: Description of the “TRIFFID” dynamic global vegetation model, Hadley Centre Technical Note, Met Office Hadley Centre, Exeter, Devon, UK, 16 pp., 2001. a

Dunstone, N. J., Smith, D. M., Booth, B. B. B., Hermanson, L., and Eade, R.: Anthropogenic aerosol forcing of Atlantic tropical storms, Nat. Geosci., 6, 534–539,, 2013. a

Eastman, R. and Wood, R.: Factors controlling low-cloud evolution over the eastern subtropical oceans: A Lagrangian perspective using the A-train satellites, J. Atmos. Sci., 73, 331–351,, 2016. a

Eastman, R. and Wood, R.: The competing effects of stability and humidity on subtropical stratocumulus entrainment and cloud evolution from a Lagrangian perspective, J. Atmos. Sci., 75, 2563–2578,, 2018. a

Eastman, R., Wood, R., and Bretherton, C. S.: Time scales of clouds and cloud-controlling variables in subtropical stratocumulus from a Lagrangian perspective, J. Atmos. Sci., 73, 3079–3091,, 2016. a

Elsaesser, G. S., O'Dell, C. W., Lebsock, M. D., Bennartz, R., Greenwald, T. J., and Wentz, F. J.: The multisensor Advanced climatology of liquid water path (MAC-LWP), J. Climate, 30, 10193–10210,, 2017. a, b, c

Feingold, G., Koren, I., Yamaguchi, T., and Kazil, J.: On the reversibility of transitions between closed and open cellular convection, Atmos. Chem. Phys., 15, 7351–7367,, 2015. a

Furtado, K., Field, P. R., Boutle, I. A., Morcrette, C. J., and Wilkinson, J. M.: A Physically Based Subgrid Parameterization for the Production and Maintenance of Mixed-Phase Clouds in a General Circulation Model, J. Atmos. Sci., 73, 279–291,, 2016. a, b

Ghan, S. J.: Technical Note: Estimating aerosol effects on cloud radiative forcing, Atmos. Chem. Phys., 13, 9971–9974,, 2013. a

Golaz, J. C., Golaz, J. C., and Levy, H.: Cloud tuning in a coupled climate model: Impact on 20th century warming, Geophys. Res. Lett., 40, 2246–2251,, 2013. a

Gordon, H., Field, P. R., Abel, S. J., Dalvi, M., Grosvenor, D. P., Hill, A. A., Johnson, B. T., Miltenberger, A. K., Yoshioka, M., and Carslaw, K. S.: Large simulated radiative effects of smoke in the south-east Atlantic, Atmos. Chem. Phys., 18, 15261–15289,, 2018. a

Grosvenor, D. P. and Wood, R.: The effect of solar zenith angle on MODIS cloud optical and microphysical retrievals within marine liquid water clouds, Atmos. Chem. Phys., 14, 7291–7321,, 2014. a

Grosvenor, D. P. and Wood, R.: Daily MODIS (MODerate Imaging Spectroradiometer) derived cloud droplet number concentration global dataset for 2003–2015, availabe at: (last access: 11 December 2020), 2018. a

Grosvenor, D. P., Field, P. R., Hill, A. A., and Shipway, B. J.: The relative importance of macrophysical and cloud albedo changes for aerosol-induced radiative effects in closed-cell stratocumulus: insight from the modelling of a case study, Atmos. Chem. Phys., 17, 5155–5183,, 2017. a, b, c

Grosvenor, D. P., Sourdeval, O., and Wood, R.: Parameterizing cloud top effective radii from satellite retrieved values, accounting for vertical photon transport: quantification and correction of the resulting bias in droplet concentration and liquid water path retrievals, Atmos. Meas. Tech., 11, 4273–4289,, 2018a. a, b, c

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., 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,, 2018b. a, b, c, d, e

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. a

Gustafson, W. I. and Yu, S.: Generalized approach for using unbiased symmetric metrics with negative values: Normalized mean bias factor and normalized mean absolute error factor, Atmos. Sci. Lett., 13, 262–267,, 2012. a, b, c, d

Han, Q., Rossow, W. B., Chou, J., and Welch, R. M.: Global variations of column droplet concentration in low-level clouds, Geophys. Res. Lett., 25, 1419–1422,, 1998. a

Hanna, E., Jones, J. M., Cappelen, J., Mernild, S. H., Wood, L., Steffen, K., and Huybrechts, P.: The influence of North Atlantic atmospheric and oceanic forcing effects on 1900–2010 Greenland summer climate and ice melt/runoff, Int. J. Climatol., 33, 862–880,, 2012. a

Hill, A. A., Feingold, G., and Jiang, H.: The Influence of Entrainment and Mixing Assumption on Aerosol-Cloud Interactions in Marine Stratocumulus, J. Atmos. Sci., 66, 1450–1464,, 2009. a, b

Hoerling, M., Hurrell, J., Eischeid, J., and Phillips, A.: Detection and Attribution of Twentieth-Century Northern and Southern African Rainfall Change, J. Climate, 19, 3989–4008,, 2006. a

Holland, D. M., Thomas, R. H., de Young, B., Ribergaard, M. H., and Lyberth, B.: Acceleration of Jakobshavn Isbrætriggered by warm subsurface ocean waters, Nat. Geosci., 1, 659–664,, 2008. a

Jiang, Y., Lu, Z., Liu, X., Qian, Y., Zhang, K., Wang, Y., and Yang, X.-Q.: Impacts of global open-fire aerosols on direct radiative, cloud and surface-albedo effects simulated with CAM5, Atmos. Chem. Phys., 16, 14805–14824,, 2016. a

Khairoutdinov, M. and Kogan, Y.: A New Cloud Physics Parameterization in a Large-Eddy Simulation Model of Marine Stratocumulus, Mon. Weather Rev., 128, 229–243,<0229:ANCPPI>2.0.CO;2, 2000. a

Knight, J. R.: A signature of persistent natural thermohaline circulation cycles in observed climate, Geophys. Res. Lett., 32, L20708,, 2005. a

Knight, J. R., Folland, C. K., and Scaife, A. A.: Climate impacts of the Atlantic Multidecadal Oscillation, Geophys. Res. Lett., 33, L17706,, 2006. a

Kuhlbrodt, T., Jones, C. G., Sellar, A., Storkey, D., Blockley, E., Stringer, M., Hill, R., Graham, T., Ridley, J., Blaker, A., Calvert, D., Copsey, D., Ellis, R., Hewitt, H., Hyder, P., Ineson, S., Mulcahy, J., Siahaan, A., and Walton, J.: The Low-Resolution Version of HadGEM3 GC3.1: Development and Evaluation for Global Climate, J. Adv. Model. Earth Syst., 10, 2865–2888,, 2018. a

Lebsock, M. and Su, H.: Application of active spaceborne remote sensing for understanding biases between passive cloud water path retrievals, J. Geophys. Res.-Atmos., 119, 8962–8979,, 2014. a, b, c

Levy, R. C., Mattoo, S., Munchak, L. A., Remer, L. A., Sayer, A. M., Patadia, F., and Hsu, N. C.: The Collection 6 MODIS aerosol products over land and ocean, Atmos. Meas. Tech., 6, 2989–3034,, 2013. a

Liu, Y., Daum, P. H., Guo, H., and Peng, Y.: Dispersion bias, dispersion effect, and the aerosol-cloud conundrum, Environ. Res. Lett., 3, 045021,, 2008. a, b, c, d

Liu, Y., Zhang, K., Qian, Y., Wang, Y., Zou, Y., Song, Y., Wan, H., Liu, X., and Yang, X.-Q.: Investigation of short-term effective radiative forcing of fire aerosols over North America using nudged hindcast ensembles, Atmos. Chem. Phys., 18, 31–47,, 2018. a

Lock, A. P.: The Numerical Representation of Entrainment in Parameterizations of Boundary Layer Turbulent Mixing, Mon. Weather Rev., 129, 1148–1163,<1148:TNROEI>2.0.CO;2, 2001. a

Lock, A. P., Brown, A. R., Bush, M. R., Martin, G. M., and Smith, R. N. B.: A New Boundary Layer Mixing Scheme. Part I: Scheme Description and Single-Column Model Tests, Mon. Wea. Rev., 128, 3187–3199,<3187:anblms>;2, 2000. a

Malavelle, F. F., Haywood, J. M., Jones, A., Gettelman, A., Clarisse, L., Bauduin, S., Allan, R. P., Karset, I. H. H., Kristjánsson, J. E., Oreopoulos, L., Cho, N., Lee, D., Bellouin, N., Boucher, O., Grosvenor, D. P., Carslaw, K. S., Dhomse, S., Mann, G. W., Schmidt, A., Coe, H., Hartley, M. E., Dalvi, M., Hill, A. A., Johnson, B. T., Johnson, C. E., Knight, J. R., O'Connor, F. M., Partridge, D. G., Stier, P., Myhre, G., Platnick, S., Stephens, G. L., Takahashi, H., and Thordarson, T.: Strong constraints on aerosol-cloud interactions from volcanic eruptions, Nature, 546, 485–491,, 2017. a, b, c, d, e, f, g

Mann, G. W., Carslaw, K. S., Spracklen, D. V., Ridley, D. A., Manktelow, P. T., Chipperfield, M. P., Pickering, S. J., and Johnson, C. E.: Description and evaluation of GLOMAP-mode: a modal global aerosol microphysics model for the UKCA composition-climate model, Geosci. Model Dev., 3, 519–551,, 2010. a, b

Mann, G. W., Carslaw, K. S., Ridley, D. A., Spracklen, D. V., Pringle, K. J., Merikanto, J., Korhonen, H., Schwarz, J. P., Lee, L. A., Manktelow, P. T., Woodhouse, M. T., Schmidt, A., Breider, T. J., Emmerson, K. M., Reddington, C. L., Chipperfield, M. P., and Pickering, S. J.: Intercomparison of modal and sectional aerosol microphysics representations within the same 3-D global chemical transport model, Atmos. Chem. Phys., 12, 4449–4476,, 2012. a, b

Martin, G. M., Johnson, D. W., and Spice, A.: The measurement and parameterization of effective radius of droplets in warm stratocumulus clouds, J. Atmos. Sci., 51, 1823–1842,<1823:TMAPOE>2.0.CO;2, 1994. a

McCarthy, G. D., Haigh, I. D., Hirschi, J. J.-M., Grist, J. P., and Smeed, D. A.: Ocean impact on decadal Atlantic climate variability revealed by sea-level observations, Nature, 521, 508–510,, 2015. a

McCoy, D. T., Field, P., Gordon, H., Elsaesser, G. S., and Grosvenor, D. P.: Untangling causality in midlatitude aerosol–cloud adjustments, Atmos. Chem. Phys., 20, 4085–4103,, 2020. a

McCoy, D. T., Bender, F. A.-M., Grosvenor, D. P., Mohrmann, J. K., Hartmann, D. L., Wood, R., and Field, P. R.: Predicting decadal trends in cloud droplet number concentration using reanalysis and satellite data, Atmos. Chem. Phys., 18, 2035–2047,, 2018. a

Menary, M. B., Hodson, D. L. R., Robson, J. I., Sutton, R. T., Wood, R. A., Menary, M. B., Hodson, D. L. R., Robson, J. I., Sutton, R. T., and Wood, R. A.: A Mechanism of Internal Decadal Atlantic Ocean Variability in a High-Resolution Coupled Climate Model, J. Climate, 28, 7764–7785,, 2015. a

Miles, N. L., Verlinde, J., and Clothiaux, E. E.: Cloud droplet size distributions in low-level stratiform clouds, J. Atmos. Sci., 57, 295–311,<0295:CDSDIL>2.0.CO;2, 2000. a

Miltenberger, A. K., Field, P. R., Hill, A. A., Rosenberg, P., Shipway, B. J., Wilkinson, J. M., Scovell, R., and Blyth, A. M.: Aerosol–cloud interactions in mixed-phase convective clouds – Part 1: Aerosol perturbations, Atmos. Chem. Phys., 18, 3119–3145,, 2018a. a

Miltenberger, A. K., Field, P. R., Hill, A. A., Shipway, B. J., and Wilkinson, J. M.: Aerosol–cloud interactions in mixed-phase convective clouds – Part 2: Meteorological ensemble, Atmos. Chem. Phys., 18, 10593–10613,, 2018b. a

Morcrette, C. J.: Improvements to a prognostic cloud scheme through changes to its cloud erosion parametrization, Atmos. Sci. Lett., 13, 95–102,, 2012. a

Mulcahy, J. P., Johnson, C., Jones, C. G., Povey, A. C., Scott, C. E., Sellar, A., Turnock, S. T., Woodhouse, M. T., Abraham, N. L., Andrews, M. B., Bellouin, N., Browse, J., Carslaw, K. S., Dalvi, M., Folberth, G. A., Glover, M., Grosvenor, D., Hardacre, C., Hill, R., Johnson, B., Jones, A., Kipling, Z., Mann, G., Mollard, J., O'Connor, F. M., Palmieri, J., Reddington, C., Rumbold, S. T., Richardson, M., Schutgens, N. A. J., Stier, P., Stringer, M., Tang, Y., Walton, J., Woodward, S., and Yool, A.: Description and evaluation of aerosol in UKESM1 and HadGEM3-GC3.1 CMIP6 historical simulations, Geosci. Model Dev. Discuss.,, in review, 2020. a, b, c, d, e, f

Mulcahy, J. P., Jones, C., Sellar, A., Johnson, B., Boutle, I. A., Jones, A., Andrews, T., Rumbold, S. T., Mollard, J., Bellouin, N., Johnson, C. E., Williams, K. D., Grosvenor, D. P., and McCoy, D. T.: Improved Aerosol Processes and Effective Radiative Forcing in HadGEM3 and UKESM1, J. Adv. Model. Earth Syst., 10, 2786–2805,, 2018. a, b, c

Mülmenstädt, J. and Feingold, G.: The Radiative Forcing of Aerosol-Cloud Interactions in Liquid Clouds: Wrestling and Embracing Uncertainty, Curr. Clim. Chang. Reports, 4, 23–40,, 2018. a

Mülmenstädt, J., Gryspeerdt, E., Salzmann, M., Ma, P.-L., Dipu, S., and Quaas, J.: Separating radiative forcing by aerosol–cloud interactions and rapid cloud adjustments in the ECHAM–HAMMOZ aerosol–climate model using the method of partial radiative perturbations, Atmos. Chem. Phys., 19, 15415–15429,, 2019. a, b, c, d

Myhre, G., Shindell, D., Bréon, F.-M., Collins, W., Fuglestvedt, J., Huang, J., Koch, D., Lamarque, J.-F., Lee, D., Mendoza, B., Nakajima, T., Robock, A., Stephens, G., Takemura, T., and Zhang, H.: Anthropogenic and Natural Radiative Forcing, in: limate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, edited by: Stocker, T. F., Qin, D., Plattner, G.-K., Tignor, M., Allen, S. K., Boschung, J., Nauels, A., Xia, Y., Bex, V., and Midgley, P. M., Cambridge University Press, Cambridge, UK and New York, USA, 659–740,, 2013. a

Nakajima, T. and King, M. D.: Determination of the optical-thickness and effective particle radius of clouds from reflected solar-radiation measurements 1. Theory, J. Atmos. Sci., 47, 1878–1893,<1878:DOTOTA>2.0.CO;2, 1990. a

O'Dell, C. W., Wentz, F. J., and Bennartz, R.: Cloud liquid water path from satellite-based passive microwave observations: A new climatology over the global oceans, J. Clim., 21, 1721–1739,, 2008. a

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, D24206,, 2011. a, b, c

Robson, J., Ortega, P., and Sutton, R.: A reversal of climatic trends in the North Atlantic since 2005, Nat. Geosci., 9, 513–517,, 2016. a

Robson, J., Sutton, R. T., Archibald, A., Cooper, F., Christensen, M., Gray, L. J., Holliday, N. P., Macintosh, C., McMillan, M., Moat, B., Russo, M., Tilling, R., Carslaw, K., Desbruyères, D., Embury, O., Feltham, D. L., Grosvenor, D. P., Josey, S., King, B., Lewis, A., McCarthy, G. D., Merchant, C., New, A. L., O'Reilly, C. H., Osprey, S. M., Read, K., Scaife, A., Shepherd, A., Sinha, B., Smeed, D., Smith, D., Ridout, A., Woollings, T., and Yang, M.: Recent multivariate changes in the North Atlantic climate system, with a focus on 2005–2016, Int. J. Climatol., 38, 5050–5076,, 2018. a

Rossow, W. B., Schiffer, R. A., Rossow, W. B., and Schiffer, R. A.: Advances in Understanding Clouds from ISCCP, B. Am. Meteorol. Soc., 80, 2261–2287,<2261:AIUCFI>2.0.CO;2, 1999. a

Salomonson, V. V., Barnes, W. L., Maymon, P. W., Montgomery, H. E., and Ostrow, H.: MODIS: advanced facility instrument for studies of the Earth as a system, IEEE Trans. Geosci. Remote Sens., 27, 145–153,, 1998. a

Seethala, C. and Horváth, Á.: Global assessment of AMSR-E and MODIS cloud liquid water path retrievals in warm oceanic clouds, J. Geophys. Res., 115, D13202,, 2010. a

Seinfeld, J. H. and Pandis, S. N.: Atmospheric chemistry and physics : from air pollution to climate change, Wiley, Hoboken, NJ, USA, 2006. a

Sellar, A. A., Jones, C. G., Mulcahy, J., Tang, Y., Yool, A., Wiltshire, A., O'Connor, F. M., Stringer, M., Hill, R., Palmieri, J., Woodward, S., Mora, L., Kuhlbrodt, T., Rumbold, S., Kelley, D. I., Ellis, R., Johnson, C. E., Walton, J., Abraham, N. L., Andrews, M. B., Andrews, T., Archibald, A. T., Berthou, S., Burke, E., Blockley, E., Carslaw, K., Dalvi, M., Edwards, J., Folberth, G. A., Gedney, N., Griffiths, P. T., Harper, A. B., Hendry, M. A., Hewitt, A. J., Johnson, B., Jones, A., Jones, C. D., Keeble, J., Liddicoat, S., Morgenstern, O., Parker, R. J., Predoi, V., Robertson, E., Siahaan, A., Smith, R. S., Swaminathan, R., Woodhouse, M. T., Zeng, G., and Zerroukat, M.: UKESM1: Description and evaluation of the UK Earth System Model, J. Adv. Model. Earth Syst., 11, 4513–4558,, 2019. a, b, c

Smith, D. M., Eade, R., Dunstone, N. J., Fereday, D., Murphy, J. M., Pohlmann, H., and Scaife, A. A.: Skilful multi-year predictions of Atlantic hurricane frequency, Nat. Geosci., 3, 846–849,, 2010. a

Stevens, B., Cotton, W. R., Feingold, G., and Moeng, C.-H.: Large-Eddy Simulations of Strongly Precipitating, Shallow, Stratocumulus-Topped Boundary Layers, J. Atmos. Sci., 55, 3616–3638,<3616:lesosp>;2, 1998. a, b

Stevens, R. G., Loewe, K., Dearden, C., Dimitrelos, A., Possner, A., Eirund, G. K., Raatikainen, T., Hill, A. A., Shipway, B. J., Wilkinson, J., Romakkaniemi, S., Tonttila, J., Laaksonen, A., Korhonen, H., Connolly, P., Lohmann, U., Hoose, C., Ekman, A. M. L., Carslaw, K. S., and Field, P. R.: A model intercomparison of CCN-limited tenuous clouds in the high Arctic, Atmos. Chem. Phys., 18, 11041–11071,, 2018. a

Sutton, R. T. and Dong, B.: Atlantic Ocean influence on a shift in European climate in the 1990s, Nat. Geosci., 5, 788–792,, 2012. a

Sutton, R. T. and Hodson, D. L. R.: Atlantic Ocean Forcing of North American and European Summer Climate, Science, 309, 115–118,, 2005. a

Taylor, K., Williamson, D., and Zwiers, F.: The sea surface temperature and sea ice concentration boundary conditions for AMIP II simulations, PCMDI Rep. 60, Progr. Clim. Model Diagnosis Intercomparison, Lawrence Livermore Natl. Lab., 24 pp., 2000. a

Twomey, S.: The Influence of Pollution on the Shortwave Albedo of Clouds, J. Atmos. Sci., 34, 1149–1152,<1149:tiopot>;2, 1977. a

Walters, D., Baran, A. J., Boutle, I., Brooks, M., Earnshaw, P., Edwards, J., Furtado, K., Hill, P., Lock, A., Manners, J., Morcrette, C., Mulcahy, J., Sanchez, C., Smith, C., Stratton, R., Tennant, W., Tomassini, L., Van Weverberg, K., Vosper, S., Willett, M., Browse, J., Bushell, A., Carslaw, K., Dalvi, M., Essery, R., Gedney, N., Hardiman, S., Johnson, B., Johnson, C., Jones, A., Jones, C., Mann, G., Milton, S., Rumbold, H., Sellar, A., Ujiie, M., Whitall, M., Williams, K., and Zerroukat, M.: The Met Office Unified Model Global Atmosphere 7.0/7.1 and JULES Global Land 7.0 configurations, Geosci. Model Dev., 12, 1909–1963,, 2019. a, b

Wentz, F. J. and Meissner, T.: AMSR Ocean Algorithm, Version 2, report number 121599A-1, Tech. rep., Remote Sensing Systems, Santa Rosa, CA, availabe at: (last access: 11 December 2020), 2000. a

West, R. E. L., Stier, P., Jones, A., Johnson, C. E., Mann, G. W., Bellouin, N., Partridge, D. G., and Kipling, Z.: The importance of vertical velocity variability for estimates of the indirect aerosol effects, Atmos. Chem. Phys., 14, 6369–6393,, 2014. a

Williams, K. D., Copsey, D., Blockley, E. W., Bodas-Salcedo, A., Calvert, D., Comer, R., Davis, P., Graham, T., Hewitt, H. T., Hill, R., Hyder, P., Ineson, S., Johns, T. C., Keen, A. B., Lee, R. W., Megann, A., Milton, S. F., Rae, J. G., Roberts, M. J., Scaife, A. A., Schiemann, R., Storkey, D., Thorpe, L., Watterson, I. G., Walters, D. N., West, A., Wood, R. A., Woollings, T., and Xavier, P. K.: The Met Office Global Coupled Model 3.0 and 3.1 (GC3.0 and GC3.1) Configurations, J. Adv. Model. Earth Syst., 10, 357–380,, 2017. a

Wilson, D. R. and Ballard, S. P.: A microphysically based precipitation scheme for the UK Meteorological Office Unified Model, Q. J. R. Meteorol. Soc., 125, 1607–1636,, 1999. a, b

Wilson, D. R., Bushell, A. C., Kerr-Munslow, A. M., Price, J. D., and Morcrette, C. J.: PC2: A prognostic cloud fraction and condensation scheme, I: Scheme description, Q. J. R. Meteorol. Soc., 134, 2093–2107,, 2008a. a

Wilson, D. R., Bushell, A. C., Kerr-Munslow, A. M., Price, J. D., Morcrette, C. J., and Bodas-Salcedo, A.: PC2: A prognostic cloud fraction and condensation scheme, II: Climate model simulations, Q. J. R. Meteorol. Soc., 134, 2109–2125,, 2008b. a

Wood, R.: Drizzle in stratiform boundary layer clouds, Part 1: Vertical and horizontal structure, J. Atmos. Sci., 62, 3011–3033,, 2005. a

Wood, R.: Stratocumulus Clouds, Mon. Weather Rev., 140, 2373–2423,, 2012. a, b, c

Wood, R., Stemmler, J. D., Rémillard, J., and Jefferson, A.: Low-CCN concentration air masses over the eastern North Atlantic: Seasonality, meteorology, and drivers, J. Geophys. Res., 122, 1203–1223,, 2017. a, b

Woodward, S.: Modeling the atmospheric life cycle and radiative impact of mineral dust in the Hadley Centre climate model, J. Geophys. Res.-Atmos., 106, 18155–18166,, 2001. a

Woollings, T., Franzke, C., Hodson, D. L. R., Dong, B., Barnes, E. A., Raible, C. C., and Pinto, J. G.: Contrasting interannual and multidecadal NAO variability, Clim. Dynam., 45, 539–556,, 2015. a

Yan, X., Zhang, R., and Knutson, T. R.: Underestimated AMOC Variability and Implications for AMV and Predictability in CMIP Models, Geophys. Res. Lett., 45, 4319–4328,, 2018. a

Yool, A., Popova, E. E., and Anderson, T. R.: MEDUSA-2.0: an intermediate complexity biogeochemical model of the marine carbon cycle for climate change and ocean acidification studies, Geosci. Model Dev., 6, 1767–1811,, 2013. a

Zhang, K., Wan, H., Liu, X., Ghan, S. J., Kooperman, G. J., Ma, P.-L., Rasch, P. J., Neubauer, D., and Lohmann, U.: Technical Note: On the use of nudging for aerosol–climate model intercomparison studies, Atmos. Chem. Phys., 14, 8631–8645,, 2014. a

Zhang, R. and Delworth, T. L.: Impact of Atlantic multidecadal oscillations on India/Sahel rainfall and Atlantic hurricanes, Geophys. Res. Lett., 33, L17712,, 2006. a

Zhang, R., Delworth, T. L., Sutton, R., Hodson, D. L. R., Dixon, K. W., Held, I. M., Kushnir, Y., Marshall, J., Ming, Y., Msadek, R., Robson, J., Rosati, A. J., Ting, M., and Vecchi, G. A.: Have Aerosols Caused the Observed Atlantic Multidecadal Variability?, J. Atmos. Sci., 70, 1135–1144,, 2013.  a

Zhang, R., Sutton, R., Danabasoglu, G., Delworth, T. L., Kim, W. M., Robson, J., and Yeager, S. G.: Comment on “The Atlantic Multidecadal Oscillation without a role for ocean circulation”, Science, 352, 1527,, 2016. a

Zhang, Z., Ackerman, A. S., Feingold, G., Platnick, S., Pincus, R., and Xue, H.: Effects of cloud horizontal inhomogeneity and drizzle on remote sensing of cloud droplet effective radius: Case studies based on large-eddy simulations, J. Geophys. Res.-Atmos., 117, 1–18,, 2012. a

Zhang, Z., Werner, F., Cho, H. M., Wind, G., Platnick, S., Ackerman, A. S., Di Girolamo, L., Marshak, A., and Meyer, K.: A framework based on 2-D Taylor expansion for quantifying the impacts of subpixel reflectance variance and covariance on cloud optical thickness and effective radius retrievals based on the bispectral method, J. Geophys. Res.-Atmos., 121, 7007–7025,, 2016. a

Zuidema, P., Westwater, E. R., Fairall, C., and Hazen, D.: Ship-based liquid water path estimates in marine stratocumulus, J. Geophys. Res., 110, D20206,, 2005. a

Short summary
Particles arising from human activity interact with clouds and affect how much of the Sun's energy is reflected away. Lack of understanding about how to represent this in models leads to large uncertainties in climate predictions. We quantify cloud responses to particles in the latest UK Met Office climate model over the North Atlantic Ocean, showing that, in contrast to suggestions elsewhere, increases in cloud coverage and thickness are important over large areas.
Final-revised paper