Ice nucleation by aerosols from anthropogenic pollution.

The formation of ice particles in the atmosphere strongly affects cloud properties and the climate. While mineral dust is known to be an effective ice nucleating particle, the role of aerosols from anthropogenic pollution in ice nucleation is still under debate. Here we probe the ice nucleation ability of different aerosol types by combining 11-year observations from multiple satellites and cloud-resolving model simulations. We find that, for strong convective systems, ice particle effective radius near cloud top decreases with increasing loading of polluted continental aerosols, because the ice formation is dominated by homogeneous freezing of cloud droplets that are smaller under more polluted conditions. In contrast, an increase in ice particle effective radius with polluted continental aerosols is found for moderate convection. Our model simulations suggest that this positive correlation is explained by enhanced heterogeneous ice nucleation and prolonged ice particle growth at larger aerosol loading, indicating that polluted continental aerosols contain a significant fraction of ice nucleating particles. Similar aerosol-ice relationships are observed for dust aerosols, further corroborating the ice nucleation ability of polluted continental aerosols. By catalyzing ice formation, aerosols from anthropogenic pollution could have profound impacts on cloud lifetime and radiative effect as well as precipitation efficiency.

The ice formation process determines cloud hydrometeor number and size, alters cloud fraction and lifetime, and subsequently affects the radiative balance[1,2]. Atmospheric ice formation also plays a pivotal role in global hydrological cycle, since most precipitation initiates via the ice-phase processes over land[3,4]. Ice nucleation occurs following two primary pathways: homogeneous freezing of supercooled cloud droplets or solution particles below about −37 °C[5], and heterogeneous nucleation triggered by ice nucleating particles (INPs)[6,7]. INPs possess surface properties lowering the energy barrier for deposition nucleation or droplet freezing[7-10]. Although representing only 1 out of 103-106 ambient particles at about −10 to −37 °C[6], INPs can exert important influence on cloud microphysical properties, such as ice particle number and size[2,6].

Previous studies on ice nucleating properties of aerosols mainly relied on experiments under controlled environmental conditions, using either laboratory-generated or field-collected aerosol samples[7,11,12]. Additionally, many studies attempted to deduce the ice nucleation ability of different aerosol constituents by examining chemical compositions of ice residuals, which are particles remaining after evaporation of field-collected ice hydrometeors[13,14]. Mineral dust and biological aerosols are commonly believed to dominate the INP budget[11,15], while several aerosol constituents of anthropogenic origin were also suggested to act as INPs, such as soot[16,17], glassy organic aerosols[7,18-20], metallic (especially lead-containing) aerosols[14,21], and solid ammonium sulfate[12]. For instance, Knopf et al.[18] indicated the potential contribution of anthropogenic organic particles to cirrus formation from measurements in the highly polluted environment of Mexico City. Cziczo et al.[14] revealed that a dominant component of ice residuals of cirrus clouds was metallic particles, next to mineral dust. Lacher et al.[22] found that injections from polluted continental air masses led to increases in INP concentrations in mixed-phase cloud conditions. However, anthropogenic pollution aerosols are chemically complex and may change over time and vary from region to region (due to photochemical processing, mixing state variations, phase state change, cloud recycling, etc.), and thus different studies showed contradictory ice nucleation activities[7,12,13,16]. To evaluate the overall ice nucleation ability of anthropogenic pollution aerosols over a large spatiotemporal domain, we adopt a top-down approach by assessing the bulk relationships between different types of aerosols and ice particle properties.

A number of space-borne sensors orbiting the Earth produce global continuous measurements of aerosol and cloud properties. However, a major difficulty in assessing ice nucleation activities from space is that those measurements only provide states of clouds resulting from numerous dynamical and microphysical processes (e.g., ice multiplication, rimming, coagulation, deposition, etc.), not just ice nucleation. Our recent study[23] showed that the response of ice particle effective radius (Rei) of cirrus clouds to aerosol loading could be either negative (“Twomey effect”) or positive (“anti-Twomey effect”), depending on meteorological conditions. The meteorological modulation of Rei-aerosol relationships directly results from the competition between homogeneous and heterogeneous ice nucleation, thus providing a new perspective to assess the role of various aerosol types in ice nucleation. In this study, we combine 11-year continuous satellite observations from the A-Train constellation and cloud-resolving simulations to scrutinize ice nucleation initiated by various aerosol types, with a focus on aerosols from anthropogenic pollution. The study region is East Asia (Fig. S1) which has a wide range of aerosol type and loading[24,25].

Results and Discussion

Satellite-derived relationships between aerosol and Rei.

We analyze changes in Rei retrieved by MODIS (Moderate Resolution Imaging Spectroradiometer) as a function of aerosol optical depth (AOD) which is a proxy of aerosol loading (Fig. 1). We focus on convective clouds with cloud-top temperature colder than −37 °C (cold-top convective clouds hereafter) and anvil cirrus generated from them. Such anvil cirrus was also classified as liquid origin cirrus by Krämer et al.[26] The temperature threshold applied here is to include the influence of homogeneous freezing. Since the retrieval of Rei by MODIS is dominated by ice particles near cloud top[27], the Rei used in this study generally corresponds to the ice particles located above the height of −37 °C isotherm. We divide all samples into three equal subsets based on cloud top height (CTH) which approximately indicates relative strength of convection[28]. We show the Rei-AOD relationships separately for the > 67th percentile group (CTH > 13.1 km, Fig. 1a) and < 33th percentile group (CTH < 11.0 km, Fig. 1b). The Rei decreases significantly with an increase in total AOD (black dashed lines in Fig. 1) for clouds with high CTH, similar to the conventional “Twomey effect” for liquid clouds[29]. In contrast, for clouds with low CTH, Rei increases remarkably with AOD at small-to-moderate AOD range (< ~0.4) and levels off at larger AOD. Similarly, we group the samples based on surface-based convective available potential energy (CAPE) which is a measure of maximal energy available for convection to consume[30]. We find a good correlation between CAPE and CTH for our samples (coefficient about 0.6). Figures 1c,d show that the Rei-AOD relationships under high and low CAPE highly resemble those for high and low CTH, respectively. Therefore, convective strength modulates the observed Rei-AOD relationships.

Figure 1.

Relationships between column AOD and Rei of cold-top convective clouds and anvil cirrus with different ranges of CTH or CAPE. (a) > 67% percentile of CTH; (b) < 33% percentile of CTH; (c) > 67% percentile of CAPE; (d) < 33% percentile of CAPE. For each CTH/CAPE group, AOD of each aerosol type is divided into four bins with increasing AOD values. The error bars denote standard error of the bin average, where σ is standard deviation and N is the sample number.

We employ satellite measurements from CALIOP (Cloud-Aerosol Lidar with Orthogonal Polarization) and AIRS (Atmospheric Infrared Sounder) to classify aerosol types. The Rei-AOD relationships for dust, polluted continental aerosol, and polluted dust (colored solid lines in Fig. 1) are similar to those for the total AOD. There are significant negative and positive correlations with relatively strong convection (implied by high CTH or CAPE) and moderate convection (low CTH or CAPE), respectively, for any of the three aerosol types. Note that polluted continental aerosols in this study do not include soil/anthropogenic dust or smoke which are classified as separate types, according to the CALIOP and AIRS retrieval algorithms[31,32]. We further separately analyze cold-top convective clouds and convection-generated anvil cirrus, and find that the opposite Rei-AOD relationships are present for both cloud types (Fig. S2). For warm-top convective clouds (cloud top temperature > −30 °C), however, Rei is insensitive to AOD of any aerosol type (Fig. S3). Here we hypothesize that aerosol effect on Rei depends on relative importance of heterogeneous and homogeneous ice nucleation which is further subject to the convective strength. To test such a hypothesis and understand underlying physical mechanisms, we perform model simulations with different perturbations of cloud condensation nuclei (CCN) and INP (see next section).

One possible cause of the observed Rei-AOD relationships is the covariations of meteorological conditions that lead to simultaneous changes in Rei and AOD. Here we use a partial correlation analysis to examine such a possibility. The partial correlation is a measure of the dependence between two variables (Rei and AOD) when the influence of possible controlling variables (meteorological parameters) is removed[33,34]. Figures 2 and S4 summarize the total correlations as well as the partial correlations with the effects of 13 meteorological parameters eliminated individually and simultaneously. Similarity between the partial and total correlations indicates that the observed Rei-AOD correlations are not significantly attributed to meteorological covariations. For each aerosol type and CTH/CAPE range, we find that the partial correlations eliminating any or all meteorological parameter(s) always have the same sign as the corresponding total correlations, and the relative differences between the two are always within 30%, indicating that the majority (≥ 70%) of the correlations is attributed to the aerosol effect. More details of our statistical method and result interpretation are provided in the Supplementary Information.

Figure 2.

Pearson’s total and partial correlations between AOD and Rei. (a) > 67% percentile of CTH; (b) < 33% percentile of CTH. The leftmost column represents the total correlation, and the other columns represent partial correlations with effects of 13 meteorological parameters eliminated all simultaneously (rightmost column) and individually (in-between columns). The magnitude of the correlation is not large, since Rei is affected by many factors other than aerosols. However, all correlations in this figure are statistically significant at the 0.01 level based on the Student’s t-test. AOD range is [0, 0.8]. The meanings of meteorological parameters are provided in Methods.

Cloud resolving simulations of the aerosol effects.

We conduct model simulations using the Weather Research and Forecasting model equipped with a spectral-bin cloud microphysics (WRF-SBM) and aerosol-aware heterogeneous ice nucleation schemes[35-37] for two cold-top convective cloud systems occurring in our satellite-analysis domain. The two cases were chosen to represent two typical convection systems with moderate and strong convective strength, respectively. For each cloud system, we perform a baseline experiment with four numerical simulations in which the number concentrations of aerosols (here confined to be hygroscopic particles with diameter > 0.1 μm that could serve as CCN under favorable supersaturation) are 100, 600, 1600, and 3200 cm−3, covering clean to polluted conditions[38-40]. In each simulation, the average INP fraction is estimated to be about 1/20000, based on heterogeneously nucleated ice number concentration predicted dynamically and the total aerosol concentration initialized in our simulations. Thus, the INP number concentration increases proportionally from about 5 to 160 L−1 when CCN increases from 100 to 3200 cm−3.

Fig. 3 summarizes the simulated Rei and fraction of heterogeneously formed ice crystals as a function of aerosol concentrations. In a rising air parcel, heterogeneous nucleation occurs ahead of homogeneous freezing because of a lower requirement on temperature/supersaturation. With strong convective strength (Fig. 3a,b), an abundant amount of water vapor is transported to the upper troposphere rapidly. Under this situation, heterogeneous nucleation consumes only a small portion of available water vapor or droplets, since INPs represent only 1 out of 20000 ambient particles. As a result, over 90% of the ice crystals are produced by homogeneous freezing of cloud droplets in all simulations (Fig. 3b). An increase in aerosol concentration leads to the formation of more and smaller cloud droplets (the “Twomey effect”[29]), and hence more and smaller ice crystals. In other words, Rei decreases with an increase in aerosol loading (Fig. 3a), consistent with the observations (Fig. 1a,c). Under moderate convective conditions (Fig. 3c,d), however, a limited amount of water vapor is available for cloud formation in the upper troposphere. The heterogeneous nucleation and subsequent ice crystal growth can efficiently compete with and even prevent homogeneous freezing, as indicated by the fact that the fraction of heterogeneously formed ice particles increases from < 10% to 50% with an increase in INPs from clean to pollution conditions (red line in Fig. 3d). Since the outburst of homogenous freezing could produce a much larger number of ice particles than its heterogeneous counterpart, an increase in INP concentration (along with an increase in aerosol concentration) would result in a net decrease in ice particle number concentration and thus an increase in Rei (red line in Fig. 3c), also consistent with the observations (Fig. 1b,d). This phenomenon, known as the “anti-Twomey effect”, has been found in previous studies[41,42]. The convective-strength modulation works not only by changing upper tropospheric water vapor amount, but also by altering updraft velocity. A stronger convection is typically characterized as larger updraft velocity, which is favourable for homogeneous freezing by reducing the growth time of heterogeneously formed ice crystals. The above mechanism also explains the insensitive response of Rei to AOD for warm-top convective clouds (cloud top temperature > −30 °C, Fig. S3). In those cases, an excessive amount of water vapor is available for ice crystal growth, because homogeneous freezing is not triggered and only limited ice particles are generated from heterogeneous nucleation.

Figure 3.

Simulated changes in Rei and fraction of heterogeneously formed ice particle number concentration with > 0.1 μm aerosols. (a, c) Rei; (b, d) fraction of heterogeneously formed ice particle number concentration; (a, b) strong convection system; (c, d) moderate convection system. In the baseline simulations (red lines), the number concentrations of INP increase proportionally with total aerosols. Black and blue lines indicate the results from two groups of sensitivity simulations in which the INP concentrations are fixed at a minimal level (0.02 L−1) and 20 L−1, respectively. The Rei are averaged over the grid points colder than −37 °C.

To further demonstrate that INPs play a key role in reproducing the observed Rei-aerosol relationships, we perform additional two groups of sensitivity simulations in which the INP concentrations are fixed at a minimal level (0.02 L−1) and 20 L−1, for the moderate convection system (black and blue lines in Fig. 3c,d). In both groups, the concentrations of aerosol with diameter > 0.1 μm vary in the same range as the baseline simulations, i.e., 100–3200 cm−3. It turns out that Rei decreases slightly with aerosol concentration in both groups of sensitivity simulations (black and blue lines in Fig. 3c), contrary to the observed positive correlation between Rei and aerosol loading (Fig. 1b,d). This means that the observed Rei-aerosol relationships (Fig. 1) would not occur if INPs were not present, or if INPs did not change with the total aerosol concentrations.

Anthropogenic pollution aerosols as INPs.

For polluted continental aerosols, the satellite data show opposite Rei-AOD relationships between strong and moderate convection (Fig. 1). Meanwhile, cloud resolving sensitivity simulations demonstrate that the opposite Rei-aerosol relationships under different convective strength can be produced only if INP concentration is roughly proportional to total aerosol concentration. Therefore, we conclude that a portion of polluted continental aerosols can serve as INPs. If INPs were absent or remained invariant along with increases in polluted continental aerosols, the observed positive correlation between Rei and polluted continental AOD under moderate convection would be reversed, as illustrated by the comparison of different groups of model simulations (Fig. 3c). Moreover, the similar observed Rei-AOD relationships between dust (well-defined INPs) and polluted continental aerosols corroborate that a portion of polluted continental aerosols possess certain ice nucleating ability. The effects of aerosols from anthropogenic pollution, serving as both INP and CCN, on microphysical properties of cold-top convective clouds are summarized in Fig. 4.

Figure 4.

Schematic of microphysical changes in cold-top convective clouds due to an increase in loading of anthropogenic pollution aerosols. Under strong convection, ice crystals above the −37 °C isotherm are primarily produced by homogeneous freezing. An increase in aerosol loading leads to the formation of more and smaller ice crystals. Under moderate convection, with an increase in INPs from clean to pollution conditions, heterogeneous nucleation gradually overtakes homogeneous freezing, resulting in a decrease in ice particle number concentration and thus an increase in Rei. Aerosol induced changes in cloud macrophysics and dynamics, such as the aerosol invigoration effect[49,50], are not illustrated.

We note that the robustness of our conclusion for anthropogenic pollution aerosols acting as INPs is subject to the uncertainty in satellite derived aerosol type. To minimize this uncertainty, we only select samples with the same aerosol type within a 50-km radius and further filter dust out of polluted continental samples using independent retrievals from AIRS. Also note that anthropogenic pollution is more severe in stagnant conditions and urban areas, while dust concentration is usually larger in windy conditions and biological aerosol is more prevalent in remote areas[43,44]. Even if dust or biological aerosols (both being well-defined INPs[11,12]) exist in some polluted continental samples occasionally, concentrations of these contaminates unlikely increase proportionally with that of polluted continental aerosols and thus still cannot explain the observed Rei-aerosol relationships. In summary, our conclusion is not prone to the imperfection of aerosol typing in satellite retrievals.

In the present study we demonstrate that aerosols from anthropogenic pollution contain a considerable fraction of INPs using a top-down approach. The results thus provide a valuable constraint on estimates of the anthropogenic INP budget in modeling studies. An unresolved question in this study is what constituent contributes a major fraction of the INPs from anthropogenic pollution. The potential candidates are discussed in Supplementary Information. In the future, we expect more bottom-up model simulations of INPs from individual anthropogenic aerosol species, and comparison of these simulations with our top-down analysis results. Furthermore, our finding has important implication for global and regional climate studies. By acting as INPs, anthropogenic aerosols could have profound impacts on cloud lifetime and radiative effect as well as precipitation efficiency[2,6,45]. To date, only few studies have considered the heterogeneous ice nucleation by certain aerosol species from anthropogenic pollution in climate models[45-47]. Incorporation of this process can change cloud glaciation rate[48] and hence result in a different anthropogenic radiative forcing from pre-industrial to the present day. It also benefits the assessment of changes in the Earth’s hydrology cycle and the distribution of water resources.

Methods

Sources and processing of satellite retrievals.

We use collocated observations of aerosol and cloud properties from CALIOP aboard CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations), and MODIS and AIRS aboard Aqua, as summarized in Table S1 in the Supplementary Information.

In this study, we are interested in cold-top convective clouds and anvil cirrus clouds generated from them. We select single-layer clouds of these two types with valid quality assurance (QA) flags, based on the CALIOP level 2 merged aerosol and cloud layer product (05kmMLay, V4.10) at a 5 km along-track resolution. The “cloud type” flag in CALIOP 05kmMLay product differentiates 7 cloud types, among which are deep convective cloud and cirrus. We define a cloud profile as cold-top convective cloud if the cloud type flag is deep convective cloud and the cloud top temperature is colder than −37°C. The deep convective clouds identified by CALIOP are high and optically thick clouds which are presumably associated with convection[51]. It is noted that, due to the limited ability of CALIOP in penetrating clouds, these clouds may not necessarily develop from near the surface, which has been used as a criterion for deep convective clouds in some other cloud typing algorithms[52]. The selection criteria of cold-top convective clouds should not affect the conclusion of this study, because the convective-strength modulation of Rei-aerosol relationships theoretically holds as long as the competition between homogeneous and heterogeneous nucleation takes place. Similarly, the warm-top convective clouds used in the analysis of Fig. S3 are selected, when the cloud type flag is deep convective cloud and the cloud top temperature is warmer than −30°C. We use a temperature threshold of > −30 °C to minimize the impact of settlement of homogeneously formed ice particles at height above the −37 °C isotherm.

The cirrus clouds identified by CALIOP include two major types with distinct formation mechanisms: convection-generated (anvil) and in-situ formed. Following Zhao et al.[23], a cirrus cloud profile is classified as convection-generated if it is physically “connected” to deep convection profiles. Two neighboring CALIOP cloud layers are considered to be “connected” if they vertically overlap and are horizontally separated by no more than one profile (i.e., distance ≤ 5 km). Zhao et al.[23] also examined the physical properties of convection-generated and in-situ formed cirrus and demonstrated the cirrus type classification to be reasonable.

Subsequently, we match retrievals from other sensors to the CALIOP 5 km profiles, and calculate a number of aerosol/cloud properties corresponding to each cloud profile. The Rei is calculated by averaging the 1 × 1 km MODIS Rei retrievals (MYD06 product, Collection 6) for which the “cloud phase” is ice and the Rei uncertainty is smaller than 100%, within a 10-km radius of a CALIOP profile. We only include in Rei calculation the MODIS pixels that vertically overlap with the CALIOP cloud layer (layer top pressure of CALIOP cloud layer−10 hPa ≤ cloud top pressure of MODIS ice pixel ≤ layer base pressure of CALIOP cloud layer), with the purpose to minimize contamination by MODIS cloud pixels that do not belong to the same cloud layer as detected by CALIOP.

For a given aerosol type/composition, AOD is roughly proportional to aerosol loading[53-56]. Therefore, we use AOD as a proxy for a loading of aerosols interacting with convective clouds and convection-generated anvil cirrus, which has been a common practice in previous satellite-based studies on aerosol-cloud interactions[23,57-62]. We use the column integrated AOD because convective clouds are developed from lower troposphere and thus can be affected by aerosols at various altitudes. Using aerosol information from specific layers may not fully capture all aerosol effects discussed in this study. Similar to Rei, AOD is estimated using the mean of all 10 × 10 km MODIS AOD retrievals (MYD04 product, Collection 6) within a 50-km radius of a CALIOP profile. As AOD retrievals by MODIS are usually missing at cloudy scenes, we do the averaging over a relatively large area to increase the number of samples with valid AOD values, following many previous studies[57,58,62]. The spatially averaged AOD should be representative of the CALIOP cloud profile considering the large spatial length scale of 40–400 km for AOD variation.

We classify the selected profiles into various aerosol types primarily based on the CALIOP 05kmMLay product, and further refine the classification using the AIRS level 1B Infrared Radiance Products (AIRIBRAD). The aerosol types distinguished by CALIOP include dust, polluted dust, clean continental, polluted continental, smoke, clean marine, and dusty marine[63]. For each CALIOP cloud profile, if all aerosol layers in all profiles within a 50-km radius possess the same aerosol type, it is defined as an aerosol environment of that particular type. This stringent selection criterion is expected to minimize the possibility that different aerosol types are mixed up in the retrievals. To ensure that the environment of polluted continental aerosols is not contaminated by dust, we eliminate a polluted continental cloud profile from analysis if the “dust_score” from the AIRS AIRIBRAD product (which is based on infrared absorption) is larger than or equal to 380 (indicating that dust probably exists[64]) at any location within a 50-km radius of that cloud profile. Only less than 1% of all profiles with polluted continental aerosols have been removed.

Finally, to evaluate the impact of meteorology on observed aerosol-cloud relationships, we obtain a series of meteorological parameters (see Table S1) from the CALIOP 05kmAPro product (V4.10) and the Final Analysis reanalysis data product of National Centers for Environmental Prediction (NCEP) with 1° × 1° and 6-h resolutions. We match the NCEP reanalysis data at 6:00 UTC, which is closest to the satellite overpassing time (5:00 to 8:00 UTC), to CALIOP cloud profiles by determining which NCEP’s 1° × 1° grid contains a CALIOP profile. The meteorological parameters used in this study include: relative humidity at 300 hPa (RH300), relative humidity at 500 hPa (RH500), relative humidity at 850 hPa (RH850), CAPE, convective inhibition (CIN), pressure vertical velocity at 500 hPa (VV500), pressure vertical velocity at 300 hPa (VV300), U-component of wind speed at 300 hPa (U300), U-component of wind speed at 1000 hPa (U1000), V-component of wind speed at 300 hPa (V300), V-component of wind speed at 1000 hPa (V1000), vertical wind shear at potential vorticity surface of 2 × 10−6 deg K m2 kg−1 s−1 (VWSH+), and vertical wind shear at potential vorticity surface of −2 × 10−6 deg K m2 kg−1 s−1 (VWSH-). We use CAPE to indicate the maximal energy for convection to consume[30,65,66]. Although mid-latitude strong convection is not necessarily characterized by high CAPE, the large-scale instability indicated by high CAPE favors the convection development at 13:30 local time when the satellites overpass. CIN is typically lower than 30 J/kg for the samples used in the current study, so it is not considered in the classification of samples in Fig. 1.

Cloud-resolving Simulations.

In this study, we conduct cloud-resolving simulations to confirm our hypothesis that the amount of INP is the key to reproducing the satellite derived Rei-aerosol relationships under different convection strength, and to illustrate the underlying mechanisms. We employ the WRF model version 3.6 with a fast SBM cloud microphysics. The simulation domain covers a 300×300 km area with a 3 km grid resolution. The fast version of SBM incorporated in the WRF model[35-37] retains the advantages of the full SBM in Khain et al.[67] and produces cloud microphysical and dynamical structure as well as precipitation similar to the full SBM[68]. SBM uses four size distribution spectra to represent CCN as well as three hydrometeors including water drops (cloud and rain), low-density ice (ice and snow), and high-density ice (graupel and hail). Each spectrum is composed of 33 mass bins and the relationship between adjacent bins is determined by the function mk = 2*mk-1 (mk is the mass of a particle in the kth bin). We use the Bigg scheme[69] to simulate the droplet homogeneous freezing in temperature below −37°C. Also, the homogeneous ice nucleation of solution particles is treated by the Liu and Penner scheme[70].

Heterogeneous ice crystal formation in our model is linked to ice-friendly aerosol (IFA) via an ice nucleation parameterization. The number concentration of IFA is a prognostic variable in the model. The INP concentration is calculated dynamically as a function of IFA (na>0.5) and temperature (T), following the DeMott et al. (2015)[71] parameterization scheme (D2015): where cf, α, β, γ, δ are constants in our model. The temperature influence on the INP is explicitly considered in our model simulations. We estimate the INP fraction by calculating the ratio of simulated heterogeneously nucleated ice number concentration to the total aerosol concentration initialized in our simulations.

In accordance with the objective of our simulations described above, we focus on general nucleating ability of aerosols (serving as either CCN or INP) in the WRF sensitivity simulations, but not explicitly considering specific aerosol type or chemical composition. Hence, the simulated Rei-aerosol relationships are independent of specific aerosol type/composition. D2015 is an updated version of DeMott et al. (2010)[6], and includes both the laboratory data from the Aerosol Interactions and Dynamics of the Atmosphere (AIDA) chamber and field campaign data. Both of them share a similar formula, which is a function of aerosol particles larger than 0.5 μm and ambient air temperature. Since our model does not explicitly consider aerosol chemical compositions, there is no fundamental differences using those two schemes, except that the updated formula has a relatively more efficient ice nucleation rate, which can realistically simulate a dust case in a recent study[35].

We simulate two cold-top convective cloud systems in real-case settings. The two systems, which occurred in March 22, 2015 and March 29, 2009, are representative of moderate and strong convection, respectively. Figure S7 shows the vertical profiles of updraft velocity, specific humidity and temperature in the simulations of two types of convection. The maximal updraft velocity is much larger in the strong convection case than in the moderate convection case. Note that we use AOD in satellite data analysis and > 0.1 μm aerosol concentrations in modelling to represent loading of aerosols interacting with clouds. Since these two variables are both linked with the total aerosol loading, the observed and modelled Rei-aerosol relationships can be compared in this study.

Code availability

The code of WRF-SBM model is available at http://www2.mmm.ucar.edu/wrf/users/download/get_source.html. The scripts used to process the satellite data can be requested from the corresponding authors.

Data availability

The satellite and meteorology data products used in this study are publicly available at the following sites.

MODIS/Aqua MYD04 and MYD06 products: https://earthdata.nasa.gov/

CALIOP/CALIPSO 05kmMLay and 05kmAPro products: https://eosweb.larc.nasa.gov/

AIRS/Aqua AIRIBRAD product: https://disc.gsfc.nasa.gov/

NCEP Final Analysis product: https://rda.ucar.edu/datasets/ds083.2/

Other data supporting the findings of this study are available within the article and its Supplementary Information.