Gaps in our understanding of ice-nucleating particle sources exposed by global simulation of the UK Earth System Model

Changes in the availability of a subset of aerosol known as ice-nucleating particles (INPs) can substantially alter cloud microphysical and radiative properties. Despite very large spatial and temporal variability in INP properties, many climate models do not currently represent the link between (i) the global distribution of aerosols and INPs and (ii) primary ice production in clouds. Here we use the UK Earth System Model to simulate the global distribution of dust, marine-sourced, and black carbon INPs suitable for immersion-mode freezing of liquid cloud droplets over an annual cycle. The model captures the overall spatial and temporal distribution of measured INP concentrations, which is strongly influenced by the world's major mineral dust source regions. A negative bias in simulated versus measured INP concentrations at higher freezing temperatures points to incorrectly defined INP properties or a missing source of INPs. We find that the ability of the model to reproduce measured INP concentrations is greatly improved by representing dust as a mixture of mineralogical and organic ice-nucleating components, as present in many soils. To improve the agreement further, we define an optimized hypothetical parameterization of dust INP activity (n_s(T)) as a function of temperature with a logarithmic slope of −0.175 K⁻¹, which is much shallower than existing parameterizations (e.g. −0.35 K⁻¹ for the K-feldspar data of Harrison et al., 2019). The results point to a globally important role for an organic component associated with mineral dust.

1 Introduction

A subset of atmospheric aerosol particles known as ice-nucleating particles (INPs) have the potential to trigger ice formation when immersed in supercooled water droplets. This process is referred to as immersion-mode heterogeneous ice nucleation and is thought to be the dominant heterogeneous freezing mechanism in mixed-phase clouds (Murray et al., 2012; Burrows et al., 2022). Immersion-mode freezing generally operates between 0 and about −37 °C (Murray et al., 2012), below which homogeneous ice nucleation dominates (Koop and Murray, 2016; Herbert et al., 2015). Such primary ice formation in supercooled clouds is followed by several microphysical processes that can substantially alter the cloud radiative properties as well as the generation of precipitation on regional to global scales (Korolev et al., 2017; Vergara-Temprado et al., 2018b; Hawker et al., 2021; Liu et al., 2012; Vergara-Temprado et al., 2017; Shi and Liu, 2019). The importance of ice nucleation was demonstrated by Vergara-Temprado et al. (2018b), who found that defining a realistic INP population over the Southern Ocean leads to an improved simulation of the radiative properties of mixed-phase clouds. Estimates of equilibrium climate sensitivity (ECS) have been shown to be dependent on the supercooled liquid cloud fraction (Tan et al., 2016) and are sensitive to the magnitude of the cloud-phase feedback (Zelinka et al., 2020; Ceppi et al., 2017). Murray et al. (2021) hypothesize a link between the availability of INPs and the magnitude of the cloud-phase feedback, thereby linking the global distribution of INP concentrations to ECS estimates. The importance of ice formation for cloud radiative properties and feedbacks provides an incentive to include INPs as an active component of weather and climate models – i.e. treating ice formation as a process that depends on the ice-nucleating properties of transported aerosol particles.

Laboratory, field, and modelling studies have identified several components of aerosol particles that make them effective immersion-mode INPs (Murray and Liu, 2022; Burrows et al., 2022). Mineral dust is a globally important INP type due to its relatively high ice-nucleating ability and abundance (Augustin-Bauditz et al., 2014; Harrison et al., 2016; Ullrich et al., 2017; Kanji et al., 2017; Vergara-Temprado et al., 2017; Price et al., 2018; Harrison et al., 2019; Chatziparaschos et al., 2023). Most mineral dust is emitted from low-latitude arid source regions such as the Sahara and the Arabian Peninsula (Kok et al., 2021), although high-latitude cold environments are also potential sources of dust and INPs (Paramonov et al., 2018; Tobo et al., 2019; Sanchez-Marroquin et al., 2020; Meinander et al., 2022; Barr et al., 2023). Most of the mineral dust freezing behaviour below −15 °C is thought to be accounted for by the potassium-rich feldspar (K-feldspar) component of dust (Atkinson et al., 2013; Harrison et al., 2016; Vergara-Temprado et al., 2017; Harrison et al., 2019; Price et al., 2018; Harrison et al., 2022). Primary marine organic aerosol (PMOA) internally mixed with sea spray aerosol (SSA) is also a globally important type of INP (Schnell and Vali, 1976; Wilson et al., 2015; McCluskey et al., 2018a; Ickes et al., 2020; Irish et al., 2019; Zhao et al., 2021; Creamean et al., 2018). This marine-sourced INP is present in magnitudes that are comparable to dust INPs in regions far from major dust sources and is thought to be particularly relevant for the Southern Ocean (Vergara-Temprado et al., 2017, 2018b; McCluskey et al., 2018b, 2019, 2023; Zhao et al., 2021). Other particle types can contribute locally or regionally to the INP population, including primary biological aerosol particles such as fungal fragments, bacteria, or pollen (Pummer et al., 2012; Haga et al., 2013; Wex et al., 2015). Organic material associated with some mineral dusts, such as agricultural soil or alluvial material from river deltas, can also enhance the ice-nucleating ability of these dusts (Tobo et al., 2014; O'Sullivan et al., 2014; Pereira Freitas et al., 2023; Barr et al., 2023; Kawai et al., 2023; Sun et al., 2024), and volcanic ash can be transiently important (Steinke et al., 2011; Mangan et al., 2017; Maters et al., 2019, 2020). Carbonaceous material from biomass burning and anthropogenic activity has also been shown to act as an INP (Schill et al., 2016; Umo et al., 2015; Vergara-Temprado et al., 2018a). Modelling studies suggest it plays a minor role in the total INP population globally (Kanji et al., 2020; Schill et al., 2020; Vergara-Temprado et al., 2018a) but may influence INP concentrations in polluted regions alongside other anthropogenic components (Tian et al., 2022).

Global and regional models have been used to estimate the concentrations and seasonal cycles of INPs relevant for mixed-phase clouds from different sources (Lohmann and Diehl, 2006; Hoose et al., 2010; Burrows et al., 2013; Vergara-Temprado et al., 2017; Hummel et al., 2018; Huang et al., 2018; Vergara-Temprado et al., 2018a; McCluskey et al., 2019; Sanchez-Marroquin et al., 2020; Shi et al., 2022; Kawai et al., 2023; Chatziparaschos et al., 2023; Zhao et al., 2021; McCluskey et al., 2023; Schill et al., 2020). The simulated INP distributions have proven useful for understanding the relative importance of different INP types and in quantifying the impacts of INPs on clouds and climate (Vergara-Temprado et al., 2017; Shi and Liu, 2019; Shi et al., 2022; Chatziparaschos et al., 2023; DeMott et al., 2010; Hawker et al., 2021; Huang et al., 2018; Hummel et al., 2018; Kawai et al., 2023). Most models simulate dust INP concentrations based on simulated dust fields (Burrows et al., 2013; McCluskey et al., 2019; Shi et al., 2022) combined with laboratory measurements of INP activity of atmospheric dust (Niemand et al., 2012). Improved laboratory understanding of the specific minerals involved in ice nucleation, such as K-feldspar (Atkinson et al., 2013; Harrison et al., 2019) or quartz (Harrison et al., 2019), has enabled global models of INP concentrations based on these specific minerals (Vergara-Temprado et al., 2017; Chatziparaschos et al., 2023). The marine INP source has been represented by coupling either the simulated PMOA concentration (Vergara-Temprado et al., 2017; Burrows et al., 2013) or SSA surface area (McCluskey et al., 2019, 2018b) to laboratory-based parameterizations that describe their ice-nucleating ability (where environmental samples were analysed in a laboratory). Similarly, the role of carbonaceous INPs from biomass burning (Schill et al., 2020) and anthropogenic activity (Vergara-Temprado et al., 2018a) has been studied using simulated black carbon aerosol (BC) coupled to laboratory-based parameterizations (e.g. Schill et al., 2016, and Vergara-Temprado et al., 2018a)

The Met Office Unified Model (UM) has previously been used at the regional scale to study ice nucleation. Vergara-Temprado et al. (2018b) focused on cyclonic systems over the Southern Ocean and found that the cloud radiative properties are highly sensitive to the abundance of INPs. However, INPs were prescribed in the regional model based on separate global simulations in an offline aerosol model (Vergara-Temprado et al., 2017). Hawker et al. (2021) also used a regional configuration of the UM and focused on tropical deep convective clouds over the Atlantic Ocean, and they found that the definition of the immersion-mode freezing parameterization strongly affected the micro- and macroscale cloud properties. In this study we present a global INP model relevant for immersion-mode freezing of liquid water droplets. Our study builds on the work of Vergara-Temprado et al. (2017), who simulated INPs in the global chemical transport model GLOMAP (Global Model of Aerosol Processes, Mann et al., 2010). Here we use the Global Atmosphere configuration GA7.1 (Walters et al., 2019) of the Hadley Centre Global Environment Model (HadGEM). GA7.1 is the atmosphere-only setup of the HadGEM climate model, which includes a version of GLOMAP as part of the UK Chemistry and Aerosol (UKCA) sub-model. Although the aerosol model in this study (GLOMAP) is similar to that used in Vergara-Temprado et al. (2017), the aerosol emissions, transport, and removal processes are very different. An expanded global dataset of INP measurements is used to evaluate the model and identify missing INP sources. Our study constitutes a first step towards the full inclusion of the ice-nucleation process in the various Met Office modelling systems for weather and climate.

The paper is prepared as follows. Section 2 introduces the climate model and the methods by which we represent the three classes of INPs (dust, marine-sourced, and BC). In Sect. 3 we establish the efficacy of our model by evaluating the simulated seasonal cycles of dust, SSA, and PMOA concentrations against ambient aerosol measurements. Section 4 first introduces the global distribution of INPs simulated by our model and then evaluates the simulated INP concentrations against a dataset of INP measurements; the INP model evaluation includes global and regional evaluations (Sect. 4.2.1–4.2.2), as well as a focus on marine-sourced INPs in the Southern Ocean (Sect. 4.2.3). Following the outcomes of the evaluation we next focus on identifying missing INPs in our model; Sect. 5 introduces a number of possible sources and determines which are most likely to be contributing to our model bias. In Sect. 6 we summarize the main findings of our study and bring the results together to inform the direction of future studies.

2 Methods

The INP model in this study uses simulated aerosol concentrations generated by the UK Earth System Model version 1 (UKESM1), a climate configuration of the UM comprised of several coupled models detailed below. INP concentrations are calculated from the simulated aerosol fields using surface-area and mass-dependent parameterizations of immersion-mode freezing according to the singular description (Vali, 2014). This approach assumes that the time dependence of the ice-nucleation process is of second order compared with the temperature dependence and that ice nucleation occurs instantaneously at certain temperature-dependent active sites distributed on the surface of dust and BC particles or within the volume of SSA. The hypothesis of active sites is supported by experimental evidence (Holden et al., 2019, 2021) and can be used to explain in a reproducible way the freezing behaviour of water induced by most materials (Vali, 2014; DeMott et al., 2018). Therefore, although time dependence of nucleation is observed in laboratory measurements (Herbert et al., 2014; Knopf et al., 2020), for the vast majority of warm- and mixed-phase clouds the dominant driver for primary ice nucleation is via the cooling associated with vertical ascent.

2.1 The UK Earth System Model

UKESM1 is an Earth system model (Sellar et al., 2020) developed by the UK Met Office and the Natural Environment Research Council (NERC). The physical core of UKESM1 is the Global Coupled climate model 3.1 (GC3.1) comprised of the Global Atmosphere 7.1 (GA7.1) and Global Land 7.0 (GL7.0) configurations that represent the atmosphere–land–ocean–sea-ice interactions (Mulcahy et al., 2018; Walters et al., 2019). We use the atmosphere-only configuration of UKESM1. For these simulations we use a standard Atmospheric Model Intercomparison Project (AMIP) setup, whereby the sea-surface temperatures and sea ice fields are prescribed according to observed 2-D time series for the year 2010. Other relevant land-surface and ocean-surface variables are prescribed using climatologies for the year 2010, including ocean-surface chlorophyll concentrations which are used in our representation of marine-sourced INPs and discussed further in Sect. 2.3. The configuration has an atmospheric resolution of 1.25°×1.875° in latitude and longitude (about 135 km) and 85 vertical levels. We perform a free-running simulation (i.e. not nudged to reanalysis meteorology) for 20 years with prescribed fields (e.g. sea-surface temperature, surface properties, non-interactive aerosol emissions) taken from the year 2010 and calculate a mean annual cycle relevant for the present-day climate. The representation of aerosol emissions relevant for this study will be presented in the following sections. Although an evaluation of UKESM1 has been carried out by Mulcahy et al. (2020), given that we are employing the atmosphere-only configuration we present additional model evaluation of aerosol species most relevant for this study in Sect. 3.

2.2 Representation of dust INPs

We use dust concentrations simulated by the Coupled Large-scale Aerosol Simulator for Studies in Climate (CLASSIC) mineral dust scheme (Woodward, 2011) within the UM. CLASSIC prognostically represents dust in six size sections spanning diameters of 0.0632 to 63.2 µm, with interactive emissions calculated every 30 min. Dust emission occurs when the near-surface dust friction velocity exceeds a certain threshold value, which is dependent on the texture of the soil as well as the moisture content of the top layer. In the present study, dust is emitted only from model grid boxes that contain a non-zero fraction of bare soil. We assume that the ice-nucleation ability of the mineral component of soil dusts is determined by the K-feldspar content of the particles. This assumption is valid for the majority of desert dusts found worldwide (Atkinson et al., 2013; Harrison et al., 2022, 2019; Price et al., 2018; Chatziparaschos et al., 2023; Reicher et al., 2018), although it underestimates the ice-nucleation ability of soil dusts that contain organic ice-nucleating material, including agricultural soil dusts (O'Sullivan et al., 2014; Tobo et al., 2014) and high-latitude soil dusts (Tobo et al., 2019; Sanchez-Marroquin et al., 2020; Shi et al., 2022; Meinander et al., 2022; Barr et al., 2023). These organic-mineral-containing soil dusts are discussed in Sect. 5.2. Unlike Vergara-Temprado et al. (2017), we do not simulate the emission of K-feldspar as a distinct tracer. Instead, we simplify the model by assuming that 5 % of the dust mass and surface area is K-feldspar, consistent with mass fractions in atmospheric mineral dust particles (Atkinson et al., 2013; Harrison et al., 2019) and surface samples (Jeong, 2024). Although this mass fraction varies with location (about 1 to 20 %; Harrison et al., 2019; Nickovic et al., 2012; Harrison et al., 2022; Adebiyi et al., 2023) and particle size (Tanguy Claquin, 1999; Perlwitz et al., 2015), the error in INP abundance caused by using a fixed fraction is likely to be within the uncertainty of most INP measurement techniques. The assumed fraction of K-feldspar is discussed further in Sect. 5.1.

Using modelled dust concentrations, the dust INP number concentration in air (N_dust-INP) at a specific temperature T in each bin i is calculated following the singular description (DeMott, 1995; Murray et al., 2012; Vali, 2008) using Eq. 1.

$$\begin{array}{}\text{(1)}& N_{{\text{dust-INP}}_i}\left(T\right)=N_{{\mathrm{dust}}_i}{f}_{{\mathrm{KF}}_i}(\mathrm{1}-e^{n_{\mathrm{s}}\left(T\right)\mathit{\pi }{d}^{\mathrm{2}}})\end{array}$$

In this equation, $N_{{\mathrm{dust}}_i}$ is the dust number concentration in each size bin (L⁻¹ of air), $f_{{\mathrm{KF}}_i}$ is the fraction of the dust mass which corresponds to K-feldspar, n_s(T) is the temperature-dependent density of ice-nucleating active sites of pure K-feldspar, and d is the mean diameter of the bin. We define n_s(T) according to the parameterization of Harrison et al. (2019). The bins are summed to obtain the total N_dust-INP (L⁻¹ of air).

2.3 Representation of marine-sourced INPs

For aerosol emitted from the ocean surface, two methods have previously been used to calculate INP activity. One is based on the mass concentration of atmospheric PMOA concentrations (Wilson et al., 2015; Zhao et al., 2021). This approach was used by Vergara-Temprado et al. (2017), obtaining good agreement with measurements of organic carbon. The reliability of the approach depends on the PMOA concentrations in the model, which are not well defined by measurements. McCluskey et al. (2019) used a simpler approach in which the ice-nucleating ability of PMOA is assumed to be proportional to the surface area of SSA, as defined by the measurements-based parameterization of McCluskey et al. (2018b). A model based on this assumption agreed well with INP measurements over the Southern Ocean and Mace Head (Ireland) (McCluskey et al., 2019). In our simulation we use the SSA-based approach to calculate marine-sourced INPs but test both approaches in Sect. 4.2.3.

For the SSA INPs we use SSA concentrations simulated by the modal version of the Global Model of Aerosol Processes (GLOMAP-mode). GLOMAP-mode, described in full by Mann et al. (2010) and Mulcahy et al. (2020), is a two-moment aerosol microphysics scheme that prognostically represents five log-normal aerosol modes ranging from particle diameters of 1 nm to 10 µm. SSA is present in two of the soluble modes (accumulation and coarse) and is internally mixed with sulfate, BC, and organic matter. The N_marine-INP concentration is calculated following the singular description (DeMott, 1995; Murray et al., 2012; Vali, 2008) according to Eq. (2).

$$\begin{array}{}\text{(2)}& N_{{\text{marine-INP}}_i}\left(T\right)=N_{{\mathrm{SSA}}_i}(\mathrm{1}-e^{n_{\mathrm{s}}\left(T\right)\mathit{\pi }{d}^{\mathrm{2}}})\end{array}$$

In order to perform these calculations, we extract the number, mass, standard deviation, and centre of the two aerosol modes that contain SSA. We then redistribute the SSA into six bins (denoted by i) from 0.1 to 10 µm. Using six bins reduces the precision in which the size distribution of the SSA can be described, but this only leads to minor effects on the marine-sourced INP concentration when compared to the uncertainties in the parameterization of the density of ice-nucleating active sites of the SSA. Using the number concentration of SSA in each bin, n_i, we calculate the INP concentration in each bin using Eq. (2), where n_s(T) is the density of ice-nucleating active sites given by McCluskey et al. (2018b). The concentration in all the bins is summed in order to calculate the total N_marine-INP.

For the PMOA-based approach we use the model emissions of PMOA associated with SSA particles (Mulcahy et al., 2020). Briefly, the fraction of marine organic content in the SSA is correlated with the presence of chlorophyll in the sea surface and negatively correlated with the wind speed at 10 m. The sea–air flux of PMOA is calculated from the flux of SSA and the fraction of marine organic material in the SSA. Adding a source of PMOA in UKESM1 substantially improved the seasonal cycle of organic aerosol mass and cloud droplet number concentrations in the Southern Ocean (Mulcahy et al., 2020).

Although PMOA is added to the SSA particles, it is transported as aerosol organic matter along with other organic components from secondary organic aerosol and primary emissions. In order to separate the PMOA from the rest of the organic carbon, we therefore ran an additional simulation where the PMOA emissions were switched off. The difference in organic aerosol mass between the two simulations therefore represents the contribution from PMOA emissions. This method was successfully applied by Mulcahy et al. (2020) using UKESM1 to compare the modelled PMOA concentrations to measurements at Amsterdam Island. However, the organic aerosol components are not always linearly additive because of processes that influence other sources and sinks of organic aerosol. This results in errors in calculated PMOA primarily in regions where PMOA makes only a small contribution to total organic aerosol (mainly in the Northern Hemisphere). Because of these issues we restrict our analysis of PMOA INPs to the Southern Ocean boundary layer, where it is an important INP.

PMOA mass is converted into primary marine organic carbon assuming a mass ratio C $/$ PMOA $=\mathrm{1.0}/\mathrm{1.9}$ (Wilson et al., 2015). The parameterization given by Wilson et al. (2015) is then used to convert the mass concentration of primary marine organic carbon into an INP concentration using Eq. (3), where [PMOC] is the mass concentration of primary marine organic carbon (in grams) and T is temperature (in kelvin).

$$\begin{array}{}\text{(3)}& N_{\text{marine-INP,PMOA}}\left(T\right)=\left[\mathrm{PMOC}\right]e^{\mathrm{11.2186}-\mathrm{0.44593}(T-\mathrm{273.15})}\end{array}$$

Although the parameterization was developed for temperatures between −6 and −27 °C, we use it in the same range as the other parameterizations (−10 to −30 °C). This method relates the bulk mass of PMOA with a number of INPs. Given that PMOA is internally mixed with sea salt in SSA, for each aerosol mode we restrict N_{marine-INP,PMOA} using the PMOA method to the number concentration of SSA.

2.4 Representation of black carbon INPs

In this model we represent INPs from biomass burning and anthropogenic activity using BC aerosol concentrations, which are simulated by GLOMAP-mode. Prescribed BC emissions in GLOMAP-mode include monthly-varying fluxes from biomass burning following van Marle et al. (2017) and annual mean fluxes from anthropogenic sources following Hoesly et al. (2018). BC is present in three soluble modes (Aitken, accumulation, and coarse) and one insoluble mode (Aitken). The N_BC-INP concentration is calculated following the singular description according to Eq. (4).

$$\begin{array}{}\text{(4)}& N_{{\text{BC-INP}}_i}\left(T\right)=N_{{\mathrm{BC}}_i}\left(\mathrm{1}-e^{n_{\mathrm{s}}\left(T\right)\mathit{\pi }{d}^{\mathrm{2}}}\right)\end{array}$$

Here $N_{{\mathrm{BC}}_i}$ is the BC number concentration in each size bin (ranging from 0.1 to 10 µm), with contributions from each aerosol size mode redistributed in the same manner as SSA in Eq. (2). We represent the ice-nucleating ability of BC aerosols using the n_s(T) parameterization presented by Schill et al. (2020). This is based on INP measurements of wildfire emissions and is consistent with diesel exhaust (Schill et al., 2016) and soot (Vergara-Temprado et al., 2018a; Ullrich et al., 2017). The n_s(T) values for wildfire emissions are slightly higher than diesel exhaust (Schill et al., 2020), so they can be seen as an upper estimate to the ice-nucleating ability of BC aerosol from all sources of simulated BC.

3 Aerosol model evaluation

In this section we compare the simulated aerosol properties relevant for immersion-mode ice nucleation with in situ measurements from ground stations. We perform this analysis for dust, SSA, and PMOA. An in-depth evaluation of aerosol in UKESM1 was performed by Mulcahy et al. (2020). Using long-term ground-based observation networks (EMEP, IMPROVE, and EANET) Mulcahy et al. (2020) demonstrated that UKESM1 reproduces the global distribution and present-day trends of aerosol from both natural and anthropogenic sources with reasonable skill. A comparison with multiple satellite products shows good agreement with the spatial distribution and seasonal cycle of aerosol optical depth, although there is low bias close to dust source regions. Checa-Garcia et al. (2021) evaluated a number of Earth system models including UKESM1 with a specific focus on dust. The authors show that UKESM1 is able to generally capture the observed deposition flux, surface concentrations, and dust loading but under-represents the dust optical depth. For our INP model we are primarily interested in representing the dust concentrations. Therefore, the lack of skill in reproducing observed dust optical depths is of secondary importance to our analysis. BC mass concentrations simulated by UKESM1 have been compared to a network of ground-based observations by Mulcahy et al. (2020). The correlation coefficient between simulated BC and measured BC is 0.45 with a slightly negative bias. This model skill is consistent with other climate models (e.g. Tegen et al., 2019, and He et al., 2015). Therefore, we expect some associated uncertainty with our simulated BC INP concentrations, which we revisit later in the paper. We extend the UKESM1 evaluations performed by Mulcahy et al. (2020) and Checa-Garcia et al. (2021) with some additional focus on the three primary sources of INP types relevant for this study: dust (Fig. 1), SSA (Fig. 2), and PMOA (Fig. 3).

Figure 1 Comparison of simulated dust concentrations with measurements. Lines correspond to the simulated monthly mean concentration at the location, and dots correspond to the monthly mean measurements from the Huneeus et al. (2011) database at 12 locations with standard deviations. The RMSE and MBE for each location are shown above each plot. The number shown in the top right of each plot is the mean number of measurements per month. Note that all stations have different y-axis scales, with Mawson, Tabuaeran, and American Samoa stations shown on a log scale.

The dataset of Huneeus et al. (2011) provides decadal-scale measurements of surface-level dust concentrations from several locations around the globe from the 1980s to the 1990s. Monthly mean measurements are shown alongside monthly mean simulated dust concentrations from the surface level of the model in Fig. 1. Overall, the model reproduces the relative magnitudes of dust concentrations and is able to capture seasonal variability. Model–measurement agreement is also good at remote locations far from dust sources, such as the Antarctic (Mawson) and the South Pacific Ocean (American Samoa and Norfolk Island). At Norfolk Island there is an inter-annually variable peak in September which is well reproduced, as are the low magnitudes (<10 µg m⁻³) at the other two stations. The greatest concentrations (20 to 100 µg m⁻³) are measured at Izana, Barbados, Bermuda, and Miami, which are all downwind from African dust source regions (Glaccum and Prospero, 1980; Prospero et al., 1981). Here, the model is able to capture the magnitudes and seasonal variability, though for winter months (November to February) the model over-predicts the dust concentrations in Izana by a factor of 4. This region, off the west coast of Africa, is subject to periodic plumes of dust that are strongly influenced by decadal-scale oscillations during the wintertime (Nakamae and Shiotani, 2013). Additionally our simulation uses some prescribed climatologies for the year 2010 so some discrepancy could be attributed to inter-annual or multi-decadal variability. A similar over-prediction (factor of 4) is observed in the tropical North Pacific island of Enewetak, influenced by the long-range transport of dust from Asia (Duce et al., 1980); similar factors could explain the discrepancies. Additional reasons for over-prediction are poorly represented dust removal processes and/or the insoluble nature of the simulated dust which would inhibit removal via cloud droplet activation. We revisit the representation and simulation of dust emissions in Sect. 5.

Figure 2 Comparison of simulated sea spray aerosol (SSA) concentrations with measurements. Lines correspond to the simulated monthly mean concentration at the location (both accumulation and coarse modes), and dots correspond to monthly mean measurements. The RMSE and MBE for each location are shown above each plot. The number shown in the top right of each plot is the mean number of measurements per month. Note the different y-axis scales used for the Mawson and Heimaey stations.

The SSA evaluation is shown in Fig. 2. Monthly mean concentrations of accumulation and coarse modes are compared with historical measurements (1980s to 1990s) from the AEROCE (Atmosphere/Ocean Chemistry Experiment) and SEAREX (Sea/Air Exchange) programmes, available from the AeroCom benchmark data website (https://aerocom.met.no/data, last access: 26 January 2024). Compared to the dust, SSA measurements display considerably less variability amongst the 12 stations, with less pronounced seasonal variation in general. The model is within the range of measurement uncertainty in most stations, although it tends to over-predict at Chatham Island and Baring Head (both situated in the region of New Zealand). Although within measurement uncertainty, the model tends to consistently under-predict the magnitude at Heimaey (North Atlantic) throughout the annual cycle by up to a factor of 4. This is important as this station has the highest concentrations of all the stations, although it also displays the greatest inter-annual variability, which may account for the model discrepancy.

Figure 3 Comparison of simulated primary marine organic aerosol (PMOA) concentrations in our simulation with measurements at Amsterdam Island. The line corresponds to the simulated monthly mean concentration at the location, and the dots correspond to monthly mean measurements. The RMSE and MBE are shown above the plot.

We compare simulated PMOA concentrations to measurements from Amsterdam Island (southern Indian Ocean) from Sciare et al. (2009) in Fig. 3. The discrepancies between the modelled and measured PMOA are consistent with Mulcahy et al. (2020), who also found that UKESM1 captures the magnitude and seasonality of the PMOA measurements in Amsterdam Island, but with the August peak concentration occurring about 2 months earlier in the model than in the measurements. This may be related to an early peak in simulated chlorophyll evident in Mulcahy et al. (2020), which is used to determine the PMOA emission flux. Nonetheless, the comparison suggests the model is able to capture the correct magnitude in PMOA, with the added caveat that there is a distinct lack of long-term PMOA measurements available to evaluate the model.

4 The INP model

4.1 INP spatial distribution

The annual mean horizontal distribution of simulated INP number concentration from the sum of marine, dust, and BC sources ($N_{\mathrm{INP}}=N_{\text{dust-INP}}+N_{\text{marine-INP}}+N_{\text{BC-INP}}$) at an altitude of 500 m is shown in Fig. 4. A reference activation temperature of −20 °C is used to calculate N_INP(T), which is the concentration that would be measured in a droplet freezing assay or an online INP counter at −20 °C, or if this air parcel were cooled to −20 °C. Dust INP concentrations (Fig. 4b) are greatest close to the major dust source regions of central and northern Africa, across much of Asia, and Australia. The outflow of these Northern Hemisphere (NH) dust sources results in elevated concentrations over the North Atlantic Ocean and North Pacific Ocean, whilst the South Pacific Ocean and Southern Ocean are influenced by outflow from Australia. There is a clear hemispheric gradient, with dust INP concentrations about 2 orders of magnitude higher in the northern high latitudes than in the south. Marine-sourced INP concentrations (Fig. 4d) are lower and less spatially variable than dust, with N_INP(−20 °C) around 10⁻³ L⁻¹ over much of the world's oceans. BC INP concentrations (Fig. 4f) are focused in regions of intense biomass burning (South America and central Africa) and regions influenced by anthropogenic pollution (Europe, North America, South Asia, and East Asia) with a maximum N_BC-INP of ∼ 0.02 L⁻¹ over East Asia. The total INP concentrations (Fig. 4a) are dominated by dust in most regions, with marine-sourced INPs making up a substantial fraction over the Southern Ocean (Fig. 4d). BC INPs contribute most strongly over the biomass burning regions, with less influence over the polluted continental regions (Fig. 4g). The presence of dust over many of the polluted regions, especially over Asia, suppresses the impact of BC INPs which reach a maximum contribution of ∼ 10 % over East Asia. This is consistent with previous global modelling studies that suggest a similarly weak role for BC INPs (Vergara-Temprado et al., 2018a; Schill et al., 2020; Kanji et al., 2020) and observational studies of anthropogenic-sourced INPs (Adams et al., 2020; Kanji et al., 2020; Zhang et al., 2022) in the immersion mode.

Figure 4 Annual mean N_INP(T) concentrations in the boundary layer at 500 m and a temperature of −20 °C for (a) total, (b) dust, (d) marine-sourced, and (f) BC N_INP(−20 °C). The marine-sourced INPs are parameterized using SSA. Panels (c), (e), and (g) show the ratio (expressed as a percentage) of N_INP(−20 °C) to total N_INP(−20 °C) for dust, marine-sourced, and BC INPs, respectively, and it is used to highlight regions that are dominated by each source of INPs.

Figure 5 Zonal mean profiles of N_INP(T_amb) across the four seasons between the surface and an altitude of 11 km. Panels (a)–(d) show the total N_INP(T_amb) with contour lines showing the seasonal mean ambient temperature at 10 °C intervals. Data for temperatures below −40 °C are muted as homogeneous freezing is likely to dominate over heterogeneous freezing. Filled contours in panels (e)–(h) show the percentage of grid boxes (over days and longitudes) where marine sources of N_INP(T_amb) are larger than the contribution from the dust sources of N_INP(T_amb).

Figure 5 shows the INP number concentration that is active at ambient temperatures, defined as N_INP(T_amb). N_INP(T_amb) is the INP concentration that we would expect to nucleate ice if activated to cloud droplets. The highest values of N_INP(T_amb) occur at high latitudes or high altitudes where temperatures are the lowest. The NH high latitudes have relatively high concentrations throughout the lowest 4 km for most of the year, caused by the long-range transport of dust from the NH low-latitude regions. Figure 5a–d show that in the Southern Hemisphere non-negligible concentrations of N_INP(T_amb) are only present within the lowest 2 km at latitudes south of 60° S, with concentrations greater than 10² L⁻¹ for most of the year. The very cold lower-tropospheric temperatures present in this region enhance the activity of the very low concentrations of aerosols there that can act as INPs. In the tropics (30° S–30° N) INPs are active at mid-tropospheric altitudes (>6 km) where aerosol concentrations tend to be lower. Here INPs may play an important role in influencing radiative properties of deep convective clouds (Hawker et al., 2021). Figure 5e–h show that below 60° S both marine and dust sources contribute to N_INP(T_amb), with marine sources more important over the Southern Ocean and dust more important over the Antarctic landmass. N_INP(T_amb) values for BC INPs do not appreciably compete (<5 %) with either marine-sourced or dust INPs and have been left out of Fig. 5. The largest contribution of marine-sourced INPs occurs during March–May (MAM), while its minimum occurs during December–February (DJF). In contrast, the contribution of marine INP sources to N_INP(T_amb) in the NH is negligible for most months of the year, only peaking at less than 5 % between October and March. Here in the NH, non-negligible N_INP(T_amb) concentrations within the lower troposphere extend as far south as 40° N for all seasons except June–August (JJA). Although the temperatures in this region are not as cold as in the Southern Hemisphere polar region, there is a greater supply of dust from the tropical and subtropical dust sources present in the NH. The development of the Arctic polar dome, evident from SON through to MAM, results in very low temperatures and high INP concentrations close to the surface. These simulations are similar to Vergara-Temprado et al. (2017), but UKESM1 simulates a much smaller contribution of marine-sourced INPs to N_INP(T_amb). This could be due to the different mechanisms used to represent marine INPs: in our study we use the SSA approach, whilst Vergara-Temprado et al. (2017) use the PMOA approach which predicts more INPs in the Southern Ocean. This is discussed further in Sect. 4.2.3.

4.2 Comparison with INP measurements

In this section we compare our INP model with a large set of ambient INP measurements made under immersion-mode conditions (at or above water saturation). We use the relevant measurements used in Vergara-Temprado et al. (2017), as well as a selection of recent measurements not included in that study (Ardon-Dryer and Levin, 2014; Price et al., 2018; McCluskey et al., 2018b; Si et al., 2018; Irish et al., 2019; Wex et al., 2019; Si et al., 2019; Jiang et al., 2020; Tobo et al., 2020; Sanchez-Marroquin et al., 2020, 2021; Tatzelt et al., 2021; Harrison et al., 2022) (see Table A1 for details; note we have not attempted to compile all available atmospheric INP data; rather we selected a subset of INP measurements with good global coverage). Figure 6 shows the location of all the measurements used in the comparison. For further analysis we split the measurements into zonal bands for the Arctic, NH mid-latitudes, NH subtropics, tropics, and the Southern Ocean. We additionally subset the Caribbean region to highlight the skill of the model in a region downwind of the Saharan dust source.

Figure 6 Map indicating the location of the INP measurements used for the evaluation of our model across the globe (a), the Arctic (b), and the Southern Ocean (c). Some markers represent multiple measurements taken from fixed stations across the measurement time period. The measurements took place across several different time periods and not only during the year of our simulation (2010). Note that the number of markers in the map is not representative of the number of INP measurements. Several measurements carried out in a fixed ground-based station are represented using a single marker, given the fact that the location does not change across the measuring period, as opposed to ship- or aircraft-based measurements. Panel (a) shows the latitudinal divisions used in this analysis, which includes the Caribbean region marked by the square box.

4.2.1 Global evaluation

The simulated N_INP and measured N_INP are compared in Fig. 7. The INP model captures the overall spatial and temporal distribution of measured INP concentrations, which shows that the model is capable of representing the sources and transport of key INP types. The majority of the measurements are within 1.5 orders of magnitude of the simulated INP concentrations, which is comparable with the spread in INP concentrations typically observed at a given location over a month (e.g. Vogel et al., 2024). However, there is some discrepancy (root mean squared error (RMSE) of 1.73 and mean bias error (MBE) of −1.31), with the model under-predicting INP concentrations, especially at lower concentrations and higher temperatures.

Figure 7b shows the ratio between the simulated N_INP and measured N_INP as a function of temperature. In general, the underestimation of the INP concentrations is greater at higher temperatures, which explains why the bias is greatest for low INP concentrations in Fig. 7a. The higher bias at low concentrations was also reported by Vergara-Temprado et al. (2017), but here we show how the bias is related to the temperature rather than the concentration. A temperature-dependent bias shows that the cause is most likely related to the assumed INP properties of the particle (i.e. INP activity versus temperature) rather than the concentration of the aerosol particles themselves, so is not related to erroneous aerosol removal rates in the model unless such errors are also very strongly temperature dependent.

In Fig. 7c we link the INP concentration bias to the n_s(T) parameterization. We express the bias as the temperature adjustment n_s(T_measurement−T_bias) that would be required to obtain the measured INP concentration. A positive bias implies that the parameterization value of n_s(T) is too low and would require a reduction in temperature of the same magnitude to higher values of n_s(T). We define this as the temperature-dependent “parameterization bias” and use the metric ΔPbias to characterize the temperature dependence of the parameterization bias for the global (or regional) dataset. ΔPbias is calculated as the difference between the mean temperature bias at temperatures higher than −17 °C and lower than −22 °C. ΔPbias is −4.3 °C for the global dataset shown in Fig. 7c and can be attributed to the bias at the higher end of the temperature spectrum, which is consistent with a missing source of INP that is active at higher temperatures. One such candidate is a contribution of ice-active organic material that is known to be associated with many soils around the world (O'Sullivan et al., 2014; Tobo et al., 2014; Steinke et al., 2016, 2020; Chen et al., 2021). Additionally, or alternatively, bias may arise from an unrepresentative parameterization of mineral dust. We revisit these topics in Sect. 5 after evaluating INP properties in different regional environments.

Figure 7 Evaluation of the INP model. (a) Simulated N_INP vs. measured N_INP. The simulated concentration shown is the monthly median for the corresponding month of the measurement. Bars correspond to the variability in the monthly medians throughout the year in the model. Diagonal lines correspond to a bias of zero, 1 order of magnitude (dashed), and 1.5 orders of magnitude, respectively. Values of RMSE and MBE in log₁₀ space are shown in the bottom-right corner. (b) Ratio between the simulated N_INP and measured N_INP versus measurement temperature. Horizontal lines correspond to a bias of zero, 1 order of magnitude, and 1.5 orders of magnitude, respectively. (c) Parameterization bias in the model for each measurement.

4.2.2 Regional evaluation

Figure 8a–f show model simulations versus measurements separated into the six regions. Here we see that the model performs better in some regions than others. Over the Southern Ocean simulated INP concentrations mostly lie within 1.5 orders of magnitude, with an overall negative bias at all concentrations. In the Arctic the model performs well at high INP concentrations but is increasingly biased low towards low concentrations. The biases are similar in the NH subtropics, which primarily includes measurements from Tokyo (Tobo et al., 2020) and Tel Aviv (Ardon-Dryer and Levin, 2014). The NH mid-latitude (western Europe and North America) measured concentrations show very little variation, with concentrations generally between 10⁻¹ and 10² L⁻¹. In contrast, the model predicts concentrations ranging from 10⁻⁶ to 10¹ L⁻¹. This represents a considerable under-prediction in our INP model's ability to represent this region. The tropics include measurements made primarily in Cabo Verde and the Amazon. The Amazon measurements (mostly at higher INP concentrations) are well represented in the model, but Cabo Verde concentrations are consistently under-predicted. Measurements from the Caribbean region are well reproduced by the INP model and largely fall along the one-to-one line. This region is strongly influenced by trans-Atlantic dust from central Africa (Harrison et al., 2022), so it is surprising that we under-predict at Cabo Verde but not downwind in the Caribbean. A potential reason for this is an unrepresentative K-feldspar fraction applied to Cabo Verde dust, which has been reported to be much greater in Cabo Verde than in the Caribbean (Glaccum and Prospero, 1980; Prospero et al., 1981). It has been reported that the K-feldspar content in Cabo Verde is around 20 % (Kandler et al., 2011), whereas in Barbados values of around 1 % have been reported (Harrison et al., 2022; Kandler et al., 2018) In addition, mineralogical maps show spatial variability in feldspar content (Nickovic et al., 2012) and that some regions (e.g. the source of Caribbean dust; the Sahel) have a feldspar content close to our assumed 5 % whilst other regions (e.g. the Sahara) have higher contents. This spatial variability in feldspar content in both source regions and of transported dust in the atmosphere could perhaps explain why we under-predict the INP concentration in Cabo Verde. This is explored further in Sect. 5.1.

Figure 8g–l show the temperature-dependent parameterization bias in each region. The regional analysis is consistent with the global dataset (Fig. 7c). At temperatures higher than −15 °C all regions display a positive bias (temperature applied to parameterizations is too high), whereas at lower temperatures the parameterization bias tends towards negative values. Each region exhibits this behaviour, but this is more apparent in some than others. In the Caribbean there is a very weak temperature-dependent parameterization bias, whereas in the Arctic and NH mid-latitudes there is a strong temperature-dependent bias. This is consistent with the regional variations in Fig. 8a–f, suggesting that in most regions the model is missing a source of INPs active at higher temperatures (above about −15 °C). The analysis suggests relatively poorer model representation at latitudes north of 23.5° N, which coincides with NH sources of pollution, including North America, Europe, and parts of Asia. In Sect. 4.1 we showed that BC contributions to the total INP concentration were secondary when compared to mineral dust. Other components of anthropogenic pollution may contribute to the INP concentration; however, studies focusing on the impact of urban pollution report weak responses (Adams et al., 2020; Chen et al., 2018; Zhang et al., 2022; Yadav et al., 2019; Bi et al., 2019). Although it is difficult to rule out an additional contribution from anthropogenic aerosols, the temperature-dependent bias is apparent in all regions except for the Caribbean, suggesting this is unlikely. The obvious INP class that is thought to be common in the terrestrial mid-latitudes and is active above $\sim -\mathrm{15}$ °C is that of biological/biogenic INPs (Kanji et al., 2017). We return to the topic of biogenic INPs in Sect. 5.2.

Figure 8 N_INP (a–f) and parameterization bias (g–l) as shown in Fig. 7a and c separated into the six regions of interest. The top row includes values of RMSE and MBE in log₁₀ space. Diagonal lines in panels (g)–(l) correspond to a bias of zero, 1 order of magnitude (dashed), and 1.5 orders of magnitude.

4.2.3 PMOA versus SSA in the Southern Ocean

In this section we focus on the effect of assuming the marine-sourced INP concentrations are related to SSA or more specifically the PMOA component of the SSA particles. We restrict the comparison to the Southern Ocean as this is the only region for which we have measurements of both PMOA and SSA (Sect. 3), and it is where our simulated marine-sourced INPs are most important (Fig. 4).

Figure 9 Comparison of simulated annual mean marine-sourced INP concentrations using the two approaches (Sect. 2.3). (a) Marine N_INP (−20 °C) concentration calculated using the PMOA approach at ground level. (b) Marine N_INP (−20 °C) concentration calculated using the SSA approach at ground level. (c) Logarithmic difference between marine N_INP (−20 °C) simulated using the PMOA and the SSA approaches (1 indicates a factor of 10 ratio). The area of the calculations has been restricted to the Southern Ocean for the reasons stated in Sect. 2.3.

Figure 9 shows the marine N_INP at −20 °C calculated using the PMOA approach (Fig. 9a) and SSA approach (Fig. 9b). Figure 9c shows the difference in the order of magnitude between the two approaches. The PMOA approach produces N_INP between 0.5 and 2 orders of magnitude higher than the SSA method. The difference between them is smaller in the northern regions of the Southern Ocean (50–60° S) as well as close to Antarctica. The spatial differences arise due to the underlying mechanisms of PMOA and SSA production, as well as removal processes in the aerosol microphysics scheme. PMOA emissions are dependent on surface concentrations of chlorophyll, whereas SSA emissions are wind speed dependent. SSA is removed from the atmosphere quicker than PMOA (Mann et al., 2010) as it is emitted into larger size modes. Differences in hygroscopicity between the two species also affect the removal rate through activation as cloud condensation nuclei. This acts to relatively enhance the transport of PMOA to the Antarctic landmass (Fig. 9c). The resulting latitudinal gradient (N–S) of the two aerosol types is more negative for SSA, driving an increasing difference in N_INP between the two approaches.

Figure 10 Comparison between simulated N_INP and measured N_INP in the Southern Ocean. The panels include (a) the measurements of McCluskey et al. (2018a) and (b–d) the measurements reported for each of the legs of the circumnavigation by Tatzelt et al. (2021) (corresponding to the Indian Ocean, Pacific Ocean, and Atlantic Ocean; see Fig. S1). The simulated N_INP, shown as lines and associated shading, is calculated for N_dust-INP+N_BC-INP (DU+BC), $N_{\text{dust-INP}}+N_{\text{BC-INP}}+N_{\text{marine-INP,SSA}}$ (DU+BC+SSA), and $N_{\text{dust-INP}}+N_{\text{BC-INP}}+N_{\text{marine-INP,PMOA}}$ (DU+BC+PMOA). Lines correspond to the median modelled concentration, while the shading corresponds to the 95 % interval over a time period and area comparable to the measurements. The measurements show the median and 25th and 75th percentiles of the data. The whiskers correspond to the 5th and 95th percentiles. Note that some data above −10 °C have been excluded from the analysis due to insufficient data at these higher temperatures.

To evaluate the two approaches we use measurements of N_INP from McCluskey et al. (2018a) and Tatzelt et al. (2021). McCluskey et al. (2018a) made INP measurements during March 2016 as part of the Clouds, Aerosols, Precipitation, Radiation, and Atmospheric Composition Over the Southern Ocean (CAPRICORN) campaign, and Tatzelt et al. (2021) made INP measurements during the Antarctic Circumnavigation Expedition (ACE). ACE took place between December 2016 and March 2017 and was split into three legs (see Fig. S1 in the Supplement): Leg 1 took place from December to January south of the Indian Ocean, Leg 2 from January to February south of the Pacific Ocean, and Leg 3 from February to March south of the Atlantic Ocean. Though INP measurements from the Southern Ocean are scarce, the spatial coverage of the combined dataset provides valuable information from across the entire Southern Ocean, albeit for a single year.

Figure 10 shows that inclusion of marine-sourced INPs, from either approach, significantly enhances N_INP when compared to dust-sourced INPs alone, especially at higher temperatures. In all cases, incorporating the marine source of INPs improves the agreement of the simulated N_INP and measured N_INP. For the McCluskey et al. (2018a) dataset (Fig. 10a), the SSA method produces a much better agreement to the measurements, whilst the PMOA-based method using the parameterization given by Wilson et al. (2015) tends to overestimate the measured concentrations. This is consistent with the results of McCluskey et al. (2019). However, for the three legs of the Tatzelt et al. (2021) dataset the approach that best matches the measurements is dependent on temperature. At lower temperatures (less than −20 °C) the measurements fall between the two approaches, whereas at higher temperatures, where most of the measurements are, the measurements are better reproduced by the PMOA-based approach. Incorrectly simulated aerosol concentrations (of SSA or PMOA) would result in a systematic bias in the simulated N_INP; therefore, the temperature dependence of the bias suggests that the cause is related to the parameterizations of marine INP activity or that the model lacks an INP aerosol species that is active at higher temperatures. We will discuss this in the following section.

The PMOA-based approach produces higher concentrations than the SSA-based approach and better reproduces the measurement dataset. However, given the limited measurements in the region, it is difficult to state which approach is ultimately more appropriate. The PMOA parameterization is based on sea-surface microlayer samples from the Arctic and the North Atlantic (Wilson et al., 2015), and inspection of that data shows that the ice-nucleating activity (per mass of organic carbon) varies by up to a factor of 100 at any one temperature, with generally lower values found for samples in the North Atlantic compared to samples from the Arctic collected in the vicinity of sea ice. Comparable measurements in the Southern Ocean are not available. The SSA-based method used in this study underestimates some of the measurements given by Tatzelt et al. (2021) but produces an INP concentration within an order of magnitude of the measurements. We conclude that the SSA method produces a valid distribution of marine N_INP over the Southern Ocean, and given it is far more simple to implement we recommend its use. It would be beneficial to repeat this analysis in other oceans and for other years, but this requires additional measurements that are currently lacking.

5 Evidence for additional INP types

The analysis in Sect. 4 suggests our INP model is unable to reproduce the temperature dependence of N_INP, with substantial bias at higher temperatures where INP concentrations are lowest. This suggests that INPs active at temperatures higher than −20 °C are either missing or currently incorrectly represented by our parameterizations. In this section we discuss potential causes of this discrepancy. Given that dust aerosol is the primary source of INPs in our model (Fig. 4), we focus on this aerosol species, though we also explore the potential for other species that are currently not simulated/represented by UKCA.

5.1 K-feldspar distribution and parameterizations

Several studies have shown that the ice-nucleating behaviour of mineral dust is primarily driven by the K-feldspar component (Atkinson et al., 2013; Harrison et al., 2016, 2019; Price et al., 2018; Harrison et al., 2022), which is the assumption we employ in this study. We have applied a constant value of 5 % K-feldspar content in our simulated dust, which is characteristic of transported dusts measured by Atkinson et al. (2013) and Harrison et al. (2019). However, we know that the K-feldspar content of soil and dust at the source is spatially variable (Nickovic et al., 2012), varying from ∼1 %–10 % (Jeong, 2024), and that larger dust particles will tend to contain relatively more K-feldspar than smaller particles (Nickovic et al., 2012). In Vergara-Temprado et al. (2017) the K-feldspar fraction of emitted dust in the model was driven by global maps of soil feldspar content and was enhanced in the larger aerosol size modes. Using this method results in a larger K-feldspar fraction close to dust source regions and a smaller fraction in remote regions. Despite the differences, the results and conclusions of that earlier study are very similar to ours, especially in the remote regions (e.g. Southern Ocean) where marine sources of INPs outcompete dust sources for much of the annual cycle. Therefore, it is reasonable to conclude that the missing high-temperature INPs in our model are not due to the spatial variability of K-feldspar that we neglect. However, the comparison between measurements in the Caribbean and the eastern tropical Atlantic indicate that there is some benefit in representing the variability in K-feldspar fraction in a similar manner to Vergara-Temprado et al. (2017).

A second source of potential model error arises from the temperature-dependent parameterization that links the dust K-feldspar content to its ice-nucleating ability. To explore this we have recalculated the global distributions of N_INP using an alternative temperature dependence of INP activity given by Atkinson et al. (2013), shown in Fig. 11a alongside the default of Harrison et al. (2019) used in our INP model. We similarly assume a K-feldspar content of 5 %. This parameterization results in higher n_s values compared to Harrison et al. (2019) across much of the temperature range but is consistent with the upper range of K-feldspar ice-nucleating activity from the compilation of literature data used by Harrison et al. (2019). The alternative K-feldspar parameterization leads to slightly better agreement between the model and the measurements (Fig. 11) in all regions except for the Caribbean, where the RMSE increases. Table 1 shows that globally the RMSE decreased from 1.73 to 1.49, ΔPbias has reduced from −4.3 to −3.0 °C, and the percentage of data with parameterization biases within a given threshold increased for all thresholds, except for the extreme cases where the parameterization bias exceeds 5 °C. This occurs because n_s is higher than the original parameterization at temperatures below −12 °C (Fig. 11a). It is not possible to say whether the Atkinson et al. (2013) parameterization is more appropriate as some of the model–measurement uncertainty biases could be due to incorrectly simulated dust transport or size distributions, or K-feldspar variability as discussed previously. Future studies would benefit from INP, dust concentrations, and K-feldspar fraction measurements made over long time periods (months to years), which would provide a means to more robustly evaluate simulated INP concentrations.

Other neglected mineralogical components may also play a role. Chatziparaschos et al. (2023) recently presented a modelling study where atmospheric dust was partitioned into components of K-feldspar and quartz. Although quartz has a lower ice-nucleating ability than K-feldspar (Harrison et al., 2019), it is more abundant than K-feldspar (Atkinson et al., 2013). Chatziparaschos et al. (2023) found that although K-feldspar dominated the overall contribution to dust-sourced INP concentrations, there were some regions and some altitudes where the quartz was more important. The greatest impact was found at low temperatures around −30 °C; therefore, the addition of this component may improve our INP model at these temperatures but is unlikely to explain the missing INPs at higher temperatures. However, the active sites on quartz samples have been shown to be unstable, with some freshly milled quartz samples becoming less active with time when exposed to both liquid water and air (Harrison et al., 2019; Kumar et al., 2019). In contrast, the active sites on K-feldspars are stable, with activity being maintained after 10 d in liquid water and a shift of only ∼1 °C after 16 months (Harrison et al., 2016). The parameterization for quartz is for freshly ground material, which might be relevant for particles that have recently undergone the saltation process, but the ageing of quartz in the atmosphere may lead to an overestimate in the contribution of quartz to the INP population.

Figure 11 Sensitivity of modelled dust INP concentrations to parameterization of ice-nucleating activity. (a) Temperature-dependent n_s curves of the parameterizations for dust INPs, where n_s for the Harrison et al. (2019) and Atkinson et al. (2013) curves is scaled by 5 % (accounting for dust K-feldspar content) to make all descriptions comparable. Panels (b), (c), and (d) show the simulated N_INP vs. measured N_INP in all regions for the Harrison et al. (2019), Atkinson et al. (2013), and O'Sullivan et al. (2014) parameterizations. Diagonal lines in panels (b)–(d) correspond to a bias of zero, 1 order of magnitude, and 1.5 orders of magnitude, respectively. The RMSE and MBE for each comparison are shown in the lower right of each panel, calculated using the log₁₀ of each dataset.

5.2 Terrestrial biological material

Previous studies have shown that terrestrial biological material contributes to the INP population (Schnell and Vali, 1976; Conen et al., 2011; Morris et al., 2014; Tobo et al., 2014; Fröhlich-Nowoisky et al., 2015; Steinke et al., 2016; Augustin-Bauditz et al., 2016; O'Sullivan et al., 2018; Steinke et al., 2020; Cornwell et al., 2023; Gong et al., 2022; Sun et al., 2024). This can be in the form of primary biological aerosol particles such as fungal fragments, bacteria, and pollen (O'Sullivan et al., 2015; Fröhlich-Nowoisky et al., 2015; Palmero et al., 2011; Pouleur et al., 1992; Morris et al., 2014; Steinke et al., 2020) or in the form of organic material, which ultimately results from biological activity, such as macromolecules or other detritus attached to agricultural soil particles (O'Sullivan et al., 2015; Conen et al., 2011; O'Sullivan et al., 2016; Hill et al., 2016; Sun et al., 2024). Recent measurements in the Arctic suggest that dust emissions from high-latitude river sediment sources have similar activity to agricultural soils and also contain an important ice-active organic component (Barr et al., 2023; Meinander et al., 2022; Tobo et al., 2019).

Figure 12 shows that the total organic carbon (TOC) content of the topsoil layer can vary by around 1 order of magnitude across relatively small spatial scales. Some dust-producing regions such as the Sahara desert have typically low TOC contents (<1 kg C m⁻²), whereas other dust emission regions such as parts of North America and central Russia have relatively high TOC contents (>10 kg C m⁻²). We might expect the contribution of ice-active organic material to dust INPs to be larger close to regions where there are periodically conditions commensurate with biological activity (i.e. liquid water), such as agricultural or seasonally fertile areas as well as river flood plains. In contrast, dust from the most arid deserts, such as much of the Sahara, has a very low organic carbon content (Minasny et al., 2017, 2014; Price et al., 2012), and measurements indicate that dusts from these arid sources are controlled by their mineral content rather than any trace of organic material (Harrison et al., 2022) (see Fig. S3).

Figure 12 Annual mean simulated dust emissions (a) and total organic carbon (TOC) content in the topsoil (0–30 cm) layer (b) taken from FAO et al. (2012).

Although some studies have attempted to include pollen grains, bacteria, and fungal spores in global models (Hoose et al., 2010; Spracklen and Heald, 2014; Cornwell et al., 2023), the lack of quantitative understanding of factors that determine fluxes of biological particles and of their ice-nucleating behaviour limits our capacity to represent all biological sources of INPs. Ice-nucleating proteins from soil fungus are known to bind to mineral dust particles and retain the ice-nucleating ability (O'Sullivan et al., 2016), whilst the membrane-bound proteins from bacteria can still nucleate ice after the cells are no longer viable and are fragmented (O'Sullivan et al., 2015; Murray et al., 2022; Bieber and Borduas-Dedekind, 2023). Additionally, pollen ice-nucleating macromolecules can be readily separated from the original pollen grains (Pummer et al., 2012; Murray et al., 2022). Hence, it has been suggested that these biological macromolecules or fragments of ice-nucleating biological entities might serve as a much more numerous group of INPs in the atmosphere than whole bacterial cells, pollen grains, or fungal spore (Rangel-Alvarado et al., 2015; O'Sullivan et al., 2015). Furthermore, the work that has been done on the global contribution of pollen, bacteria, and fungal spores to INP concentrations indicates that these species contribute very few INPs because particle concentrations in the atmosphere are much lower than concentrations of dust or SSA (Spracklen and Heald, 2014; Hoose et al., 2010).

Figure 13 Evaluation of the INP model using the fertile soil parameterization of O'Sullivan et al. (2014) to represent dust INPs. (a–f) Regional comparison between the simulated and measured concentrations and (g–l) parameterization bias as a function of measurement temperature in each region. Simulated INP concentrations shown in panels (a)–(f) are monthly medians for the corresponding month of the measurement. Diagonal lines in panels (a)–(f) correspond to a bias of zero, 1 order of magnitude (dashed), and 1.5 orders of magnitude.

Here we explore the potential effect of organic material attached to dust particles. We refer to soil dusts that contain organic material that results from biological activity as fertile soils, and we use desert soils to refer to the soils from arid regions containing very little organic matter. Fertile soils will likely produce mineral dust aerosol with more organic material associated with it than desert soils. In reality the biological content will vary between regions, e.g. Fig. 12b. Even soils in some semi-arid regions can contain organic matter and pose as a viable INP source (Schnell, 1974a, b), though airborne dust from northern Africa was found to contain no observable biological ice-nucleating components (Harrison et al., 2022). Although the global database of TOC in Fig. 12b may suggest a potential link between dust emissions and biological content, at present there are insufficient INP measurements to define the variability of activity we might expect due to the biological activity of fertile and desert soils. Therefore, as a simple first step and a reasonable approximation, we represent the ice-nucleating ability of all mineral dust as fertile soil according to the parameterization of O'Sullivan et al. (2014). The soils from the UK that O'Sullivan et al. (2014) analysed are consistent with the ice-nucleating activity of Argentinian soil dust (DeMott et al., 2018), Wyoming soil dust (Tobo et al., 2014), and high-latitude proglacial dust from the Copper River valley (Barr et al., 2023). A comparison of n_s(T) (Fig. 11a) shows that the fertile soil parameterization is 2 orders of magnitude more active than that of K-feldspar (desert soils) at temperatures higher than −12 °C but is within a factor of 10 at temperatures below −20 °C. By representing the ice-nucleating ability of all mineral dust with the ice-nucleating ability of fertile soil we expect that we will over-estimate the INP concentration in some locations (such as the Caribbean, where dust was found to contain no appreciable organic ice-nucleating material), but this test will tell us if the organic material associated with fertile soils is a potential explanation for low-biased INP concentrations elsewhere in our model.

Table 1 Table showing (top half) the RMSE and MBE between measurements and simulations using different dust INP parameterizations from the literature (Harrison et al., 2019; Atkinson et al., 2013; O'Sullivan et al., 2014) and an optimized parameterization; (middle) temperature dependence of the bias expressed using the metric ΔPbias; and (bottom half) percentage of data within parameterization bias thresholds.

Figure 11d compares the model and measurements assuming a fertile soil parameterization for all dust. A regional comparison of N_INP and parameterization bias is shown in Fig. 13. A fertile soil assumption considerably improves model–measurement agreement and removes a substantial fraction of the temperature-dependent bias in N_INP. The ΔPbias has decreased from −4.3 to −3.2 °C, and global RMSE has reduced from 1.73 to 1.01 (MBE is now 0.0), contributed by improvements in all regions except for the Caribbean (Table 1). The greatest improvement occurs in the NH subtropics and mid-latitudes, where the RMSE falls from 2.0 and 2.8 to 0.9 and 1.4. In contrast, the RMSE in the Caribbean has increased from 0.5 to 1.6. We also see an improved representation in the Arctic, where RMSE has decreased from 1.8 to 1.2. This is consistent with the simulations of Kawai et al. (2023), who found that summertime (JJA) Arctic INP measurements were best reproduced using a parameterization (Tobo et al., 2019) that accounted for measured organic content associated with locally sourced surface dust.

These fertile soil simulations suggest that the organic component in fertile soil makes an important contribution to INP activity at temperatures higher than −20 °C where INP concentrations are usually low. The results also suggest that INPs around the world could be explained by a mixture of more- or less-fertile dust particles sourced from different regions. In the Caribbean, a fully fertile soil parameterization consistently overestimates measured INP concentrations (Fig. 13f), while the desert soil assumption (Table 1) produces realistic INP concentrations. This is a region heavily influenced by the major dust source regions of central and northern Africa where dust is mostly abiotic, as demonstrated by the heat tests performed on Caribbean samples (Harrison et al., 2022) and the low TOC values at the surface (<1 kg C m⁻²; Fig. 12b). Hence, the fact that INP concentrations in the Caribbean are better represented by the desert soil parameterization is consistent with the measurements and known sources of dust. In the Arctic, the fertile soil parameterization produces INP concentrations with mostly zero bias relative to the measurements (Fig. 13h) at temperatures higher than −15 °C, but the model is generally negatively biased below −20 °C. This suggests that the Arctic is influenced by dust sourced from both fertile and desert soils. There is good evidence that dust from the low latitudes (e.g. African continent) is transported to the high latitudes alongside dust from local sources (Shi et al., 2022), the latter of which is thought to contain a high fraction of biological material (Tobo et al., 2019; Kawai et al., 2023). Therefore, we hypothesize that representing the spatial variability of K-feldspar in desert soils and organic material in fertile soils will improve the ability of our model to simulate measured INP concentrations.

5.3 High-latitude dust

A third source of uncertainty we explore is the regional pattern of dust emissions. Much of the focus of global modelling centres is to better represent the size distributions and fluxes of dust from the major dust sources, such as those in Africa, the Middle East, and Asia. Historically, less attention is given to the weaker dust source regions, which include the high latitudes. However, these regions have been identified as substantial sources of local dust (Bullard et al., 2016; Arnalds et al., 2016; Dagsson-Waldhauserova et al., 2019; Shi et al., 2022; Meinander et al., 2022), and measurements of these dust particles show high INP activity (Sanchez-Marroquin et al., 2020; Barr et al., 2023; Tobo et al., 2019; Xi et al., 2022).

Simulating high-latitude dust (HLD) is challenging given that its emission mostly occurs in glaciofluvial valleys and is strongly affected by small-scale katabatic winds (Bullard et al., 2016; Prospero et al., 2012). The driving processes tend to occur at around 10 km scales, which are much smaller than the resolution of global models and thus rely on sub-grid-scale parameterizations. Despite the challenges, some studies have simulated HLD on global and regional scales (Leadbetter et al., 2012; Groot Zwaaftink et al., 2016, 2017; Tobo et al., 2019; Hamilton et al., 2019; Mingari et al., 2020; Sanchez-Marroquin et al., 2020; Shi et al., 2022). To explore the potential impact of HLD we use UM simulations that were performed by Walters et al. (2019) using the atmosphere model identical to ours, but with additional coupling of a land-surface model. In the model description paper the authors introduce the Met Office UM Global Atmosphere and Global Land configuration (GA7.1 and GL7.0). As part of this study the authors ran simulations of the atmosphere coupled to the Joint UK Land Environment Simulator (JULES) at different resolutions including N96 (135 km resolution; as used in our INP model), N216 (60 km resolution), and N512 (25 km resolution). Although these are all unable to explicitly resolve turbulent motion, we can use them to gain a perspective on the resolution dependence of dust emissions in our model. For the following analysis we will focus on the Arctic.

Figure 14 Time series (January 2006–December 2008) of zonal mean dust emission fluxes in the NH from the high-resolution configuration (25 km) of UKESM simulations presented by Walters et al. (2019). Each coloured line shows the mean for a 10° latitude band from the Equator to 80° N. Latitudes south of 60° N are shown as solid lines, and latitudes north (the Arctic) are shown as dashed lines.

Figure 15 Dust emissions and associated surface properties during July 2008 in the UM GA7.1 and GL7.0 simulations by Walters et al. (2019) at three resolutions over the Arctic region >60° N. Monthly mean values of (a–c) dust emission fluxes, (d–f) soil moisture content, and (g–i) dust friction velocity. The configurations (N96, N216, and N512) correspond to horizontal resolutions of approximately 135, 60, and 25 km.

Figure 14 shows NH zonal mean dust emissions fluxes in the N512 resolution simulations. Peak emissions for latitudes north of 30° N occur between April and October. In the Arctic region (north of 60° N) simulated dust emissions are negligible (below 1 mg m⁻² per month) outside of these months and peak in July and August. Measurements from HLD sources (Bullard et al., 2016; Bullard and Mockford, 2018; Prospero et al., 2012) are consistent with the time period for simulated Arctic emissions (April–October) but suggest the month of peak emissions may occur earlier or later than July. This suggests even at the finest model resolution tested here UKESM may be unable to represent HLD emissions of dust from sub-grid-scale features such as valleys. Given that simulated emissions peak in July, we will focus on this month to investigate the role of model resolution dependence.

The July mean dust emission flux north of 60° N is shown for the three resolutions in Fig. 15a–c. We also show two of the primary drivers (see Sect. 4) of dust emissions in the model (Woodward et al., 2022): soil moisture content (Fig. 15d–f, shown for the top soil layer) and the dust friction velocity (Fig. 15g–i). At the coarsest resolution (135 km) dust emissions occur predominantly over North America and central Russia. As the resolution increases there are more numerous source regions, often with greater magnitudes. These new sources include coastal regions surrounding Greenland, Iceland, and the Arctic Archipelago, which likely correspond to small-scale glaciofluvial sites. Much of northern Scandinavia becomes a widespread source of dust in the finest resolution (25 km). Table 2 provides information on the July mean dust emissions and variability for the Arctic region (north of 60° N). As resolution increases, so does the annual mean flux and variability. Both the magnitude and variability (standard deviation) of the emissions increase 3-fold from the coarsest to finest resolution. This demonstrates that the overall emissions of Arctic HLD and the likelihood of extreme cases (in the form of relatively enhanced plumes of dust) are sensitive to the model resolution. This result may be specific to the model we use and its representation of dust emission processes. As an example, Shi et al. (2022) used a global model with 100 km horizontal resolution but were capable of simulating Arctic HLD emissions using a physically based dust emission model following Kok et al. (2014). Emissions in the Kok et al. (2014) dust model are dependent on wind friction velocity and the soil's threshold friction velocity, as in our model, but also include a term that accounts for increased emissions through saltation bombardment as the soil erodibility increases. Additionally, the model is primarily dependent on variables that are properties of the large-scale environment (soil and horizontal winds), which may help to suppress resolution dependence.

Table 2 Table showing the annual mean dust emissions flux and standard deviation from inside the Arctic (north of 60° N) for three different resolution configurations of the UM–UKCA coupled to JULES. The figure in the parentheses is the typical horizontal grid scale.

We explore the drivers of resolution dependence by focusing on the soil moisture content (Fig. 15d–f) and dust friction velocity (Fig. 15g–i). The soil moisture content is one of the key factors for dust uplift potential, with drier soil more readily emitted. Higher friction velocity also promotes more efficient dust uplift. The spatial distribution of soil moisture is generally similar across the scales, but the higher-resolution simulations are increasingly capable of resolving small-scale features. These include glacial valleys (visible throughout Russia and northern Asia) or sharp gradients that are smoothed out in the coarser scale (such as around the coastlines of Iceland and Greenland). Some regions have also become noticeably wetter or drier as resolution is increased, including the Alaskan North Slope and Scandinavia. The dust friction velocity displays similar sensitivity. The finer resolution permits sharper gradients and therefore more numerous extreme cases of high friction velocities. There is a particular sensitivity over the Arctic Archipelago, which drives some of the enhanced dust emissions. The reason behind these changes is beyond the scope of this paper but provides additional evidence that model resolution is important for representing Arctic HLD emissions. We hypothesize that an increase in model resolution beyond that investigated here would be required in order to reduce the discrepancy between the annual cycle of HLD emissions simulated by UKESM1 (Fig. 14) and as suggested by measurements (Prospero et al., 2012; Bullard et al., 2016; Bullard and Mockford, 2018). However, measurements of HLD are scarce, which makes evaluation difficult and provides a degree of uncertainty as to how much HLD our models should in fact simulate. This uncertainty can be observed in Shi et al. (2022), who find that although their model simulated a peak in summertime Arctic HLD at the Alert station (82.5° N, 20.3° W), this was not replicated in the measurements (Sirois and Barrie, 1999), which displayed either no seasonal cycle or a summertime minimum.

This section demonstrates that Arctic HLD emissions are sensitive to the model resolution. This may also hold true in the Southern Hemisphere, especially around the coastlines of Antarctica which are known to be potential sites for dust emissions (Asmi et al., 2018; Bory et al., 2010; Meinander et al., 2022; Kavan et al., 2017, 2018). Using a higher resolution in our simulations would likely increase the INP concentrations in our model, but this would occur primarily in the high latitudes and at all freezing temperatures. Therefore, weak HLD emissions do not explain the missing high-temperature INPs in our simulations. However, as the soil properties play a primary role in driving dust emissions, we recognize that using an interactive land-surface model would be a beneficial addition to our study and therefore recommend this for future studies focusing on HLD. Finally, there is a clear need for additional information on dust emissions from high-latitude environments, especially given that future trends suggest a shift towards more summertime melting and the appearance of either seasonal or permanent regions of bare soil that can produce dust (Meinander et al., 2022).

5.4 Optimized dust representation

Replacing the K-feldspar dust parameterization with a fertile soil assumption reduces the RMSE and some of the temperature-dependent ΔPbias (Table 1), but there remains a considerable parameterization bias towards the higher end of the temperature spectrum (Fig. 16). This suggests the INP spectrum has an even weaker temperature dependence than the fertile soils description from O'Sullivan et al. (2014). An alternative parameterization for dust INPs was produced to give the lowest combination of RMSE and ΔPbias when evaluated against the measurements. Figure 11a shows the optimized parameterization, which is described by $n_{\mathrm{s}}\left(T\right)=\exp (-\mathrm{0.4}T+\mathrm{110.777})$, where T is in kelvin, and applicable for temperatures between −10 and −30 °C. Figure 16 shows that the ΔPbias has been substantially reduced from −3.2 °C (for the O'Sullivan et al., 2014, parameterization) to −0.6 °C, and the RMSE has reduced from 1.01 to 0.71. The optimized parameterization of n_s(T) has a shallower temperature dependence than the other descriptions, with a value of −0.175 °C⁻¹ compared to −0.35 and −0.28 °C⁻¹ for the Harrison et al. (2019) and O'Sullivan et al. (2015) parameterizations. This temperature dependence is consistent with studies that have measured the ice-nucleating activity of soil and dust samples but shallower than most known INP species. Soil samples taken from four different regions of the world in Steinke et al. (2016) exhibit an n_s(T) gradient of −0.15 °C⁻¹, and dust samples from East Asian dust events in Chen et al. (2021) exhibit an n_s(T) gradient around −0.19 °C⁻¹. This suggests the material causing the shallow gradient is associated with the dust aerosol but requires further work to establish whether this is caused by organic material or another INP species. Figure 16 shows that the median parameterization bias is now close to 0 °C for measurement temperatures higher than −26 °C, but there remains a negative bias below this temperature. This supports our hypothesis that the mineral component, which exhibits a steeper temperature dependence (Fig. 11a), is important at lower temperatures, whereas the organic component (or material causing the shallow gradient) is important at higher temperatures. We note that the optimized fit is a simple exponential description that, if extrapolated to 0 °C, would be at odds with classical nucleation theory (Kashchiev, 2000). Measurements, such as those of O'Sullivan et al. (2014), demonstrate that the ice-nucleating activity of INP types often changes with temperature and rapidly decreases as the temperature approaches 0 °C. We therefore anticipate the gradient of the optimized representation, and associated ice-nucleating material, to not stay constant across the temperature spectrum.

Figure 16 Evaluation of the INP model using four representations of simulated dust INPs including three parameterizations of n_s(T) from literature (Harrison et al., 2019; Atkinson et al., 2013; O'Sullivan et al., 2014) and the optimized parameterization. (a–d) Comparison between simulated and measured INP concentrations with the RMSE for each dataset shown in the upper-left corner. (e–h) Distribution of parameterization bias (in the form of boxplots) as a function of temperature for each representation. The ΔPbias value (measure of temperature dependence) for each representation is shown in the upper-left corner of panels (e)–(h).

6 Conclusions

In this study we present the first global simulations of ice-nucleating particles (INPs) relevant for immersion-mode freezing of cloud droplets using the atmosphere configuration of the UK Earth System Model (UKESM). The model is driven by three primary sources of INPs: mineral dust aerosol assuming that its ice-nucleating ability is determined by the K-feldspar content (Harrison et al., 2019), marine sources associated with sea spray aerosol (SSA) (McCluskey et al., 2018a) or the primary marine organic aerosol (PMOA) component of the SSA (Wilson et al., 2015), and black carbon (BC) sources from natural and anthropogenic processes. Multi-decadal simulations were run in order to quantify a mean annual cycle of INPs relevant for mixed-phase clouds.

The spatial distribution of INP concentrations in UKESM is dominated by the major dust source regions, whereas marine-sourced INPs are confined primarily to the oceanic regions and have less spatial variability. BC sources influence the INP concentration close to biomass burning regions (South America, central Africa, and Indonesia) and polluted regions (North America, Europe, and Asia), but concentrations tend to be lower than dust sources. In the Northern Hemisphere (NH), dust is the dominant INP in all regions, with marine-sourced INPs playing a negligible role even over the North Atlantic Ocean, North Pacific Ocean, and Arctic Ocean. The dominance of dust in the Arctic region, albeit with a local supply of marine-sourced INPs, highlights the important role of long-range dust transport from the low-latitude dust source regions to the NH high latitudes. Dust uplifted to high altitudes contributes to the availability of high-altitude INPs that may influence high-cloud radiative properties and humidification of the stratosphere (Hawker et al., 2021; Jensen et al., 2013). In contrast to the NH, marine-sourced INPs are the dominant source of INPs over large parts of the Southern Hemisphere (SH), particularly over the Southern Ocean, where only 10 %–20 % of annual mean INP concentrations are dust. Over the Antarctic landmass, dust supplies about 60 %–80 % of INPs, which is caused by transport of dust from non-local sources (given that we have no local sources simulated in the model). This is compatible with the trajectory analysis from Neff and Bertler (2015), who found that Patagonia and New Zealand were important suppliers of dust to Antarctica, and also with McConnell et al. (2007), who have associated the desertification of Patagonia with a doubling of dust deposition to Antarctica in the twentieth century. However, other studies have identified local dust sources (Asmi et al., 2018; Bory et al., 2010) that also may play a role, especially during Austral summer melt periods.

A comparison of simulated N_INP against a large global dataset of measurements demonstrates that our INP model is able to represent the general spatial and temporal distribution of INP, which varies by roughly 4 orders of magnitude (calculated at an altitude of 500 m and −20 °C). However, the model tends to underestimate measured INP concentrations particularly when concentrations are low, which we show is caused by under-prediction at higher temperatures. As an average across all the measurements, the bias is equivalent to a 10 °C shift in the INP spectrum. Such a temperature-dependent bias occurs in all regions except for the Caribbean, which is downwind of the central Africa mineral dust sources. The tropics, NH subtropics, and NH mid-latitudes have a particularly strong temperature-dependent parameterization bias as well as an under-prediction of most measurements. The temperature dependence of the bias suggests that it cannot be due to incorrectly simulated dust transport or dust size distributions, which would cause a consistent bias across the whole temperature range.

Black carbon INPs are unlikely to be globally important but may influence regional concentrations. Over East Asia BC INPs from bio-fuel and fossil-fuel emissions represent a maximum of ∼10 % of total INPs measurable at the surface, whereas in South Asia, Europe, and North America the contribution is <5 %. Away from these regions the contribution is negligible. This is consistent with other global modelling studies (Vergara-Temprado et al., 2018a; Schill et al., 2020; Kanji et al., 2020) but suggests there may be a role for BC INPs on local scales and in synoptic conditions with low dust concentrations or in remote marine environments (Thomson et al., 2018). This is also consistent with some observational studies that suggest non-volatile organics contribute to the INP population (Tian et al., 2022; Zhao et al., 2019) but is inconsistent with others (Adams et al., 2020; Zhang et al., 2022; Yadav et al., 2019; Bi et al., 2019). The community would benefit from more compositional information on anthropogenic aerosols and the ice-nucleating activity of each component. In this study we only consider BC aerosols, yet other components, such as anthropogenic-sourced organic aerosol (Tian et al., 2022) or traffic-influenced dust emissions (Chen et al., 2024), may also act as locally important sources of INPs in polluted regions.

The INP representation and temperature-dependent model–measurement bias is improved when we assume that organic material is attached to the dust particles. The default ice-nucleation parameterization based on the K-feldspar content of dust (Harrison et al., 2019) was replaced with a parameterization based on agricultural soil samples (O'Sullivan et al., 2014) that accounts for an additional organic component. We refer to this class as fertile soil, and we use desert soil to refer to soils containing very little organic matter whose ice-nucleating ability is dependent on the mineralogical component. The fertile soil parameterization enhanced the ice-nucleating activity of dust particularly at high temperatures (greater than −15 °C), reducing some of the temperature-dependent bias (ΔPbias) from −4.3 to −3.2 °C and decreasing the RMSE of the global dataset from 1.73 to 1.01. Most improvement was found in the NH mid-latitudes and NH subtropics where the model was originally poorest. However, agreement with the Caribbean measurements is worsened and is instead best represented assuming no organic component. This demonstrates that not all dust sources can be assumed to be represented as fertile soils. We hypothesize that there is likely spatial variability in the ice-active organic component of fertile soil, which may be related to factors such as rainfall, soil organic carbon, or soil biodiversity. We aim to explore this in future studies.

An optimized parameterization of n_s(T) for dust INPs suggests the best representation exhibits a shallower temperature-dependent parameterization than current descriptions. The optimized parameterization substantially reduced ΔPbias from −3.2 °C for the O'Sullivan et al. (2014) fertile soils parameterization to −0.6 °C whilst reducing the RMSE for INP concentrations from 1.01 to 0.71. This was achieved primarily due to the weaker temperature dependence of n_s(T), which has a logarithmic slope of −0.175 °C, compared to −0.35 and −0.28 °C for the Harrison et al. (2019) and O'Sullivan et al. (2014) parameterizations. Though such gradients have been observed in laboratory studies (Steinke et al., 2016; Chen et al., 2021), further work is required to identify the aerosol material causing this and whether it is associated with the dust or is a distinct INP species that needs to be represented.

Arctic INPs are particularly important for mixed-phase clouds in cold air outbreaks, and it has been hypothesized that local sources may be important (Murray et al., 2021). We therefore investigated the model's ability to represent such sources. A comparison of simulated Arctic dust emissions using three model configurations with resolutions of 135 (this study), 60, and 25 km demonstrated that summertime high-latitude dust (HLD) emissions and their variability are considerably dependent on model resolution. Mean emissions and their variability increase 3-fold from the coarsest to the finest resolution. This suggests our INP model is missing some sources of HLD and may contribute to the misrepresentation of Arctic INP measurements. We recommend an exploration of alternative methods to reproduce HLD emissions in UKESM alongside much-needed measurements from around the Arctic in order to robustly evaluate and ground the models.

Over the Southern Ocean, we found that the simulated INP concentration is sensitive to the method used to represent the ice-nucleating ability of organics-enriched SSA. A parameterization using the PMOA content of the SSA (Vergara-Temprado et al., 2017; Burrows et al., 2013) produces more INPs than a simpler parameterization based on the SSA surface area (McCluskey et al., 2018b). The difference between the two approaches increases with proximity to the Antarctic landmass, which can be linked to differences in the emission processes of PMOA and SSA. A comparison of the two methods was made with INP measurements from the Southern Ocean. At the lowest temperatures ($\sim -\mathrm{25}$ °C) both methods produced values of N_INP consistent with the measurements, whereas at higher temperatures ($\sim -\mathrm{15}$ °C) the PMOA method is more representative. The SSA method is technically easier to represent in the model as PMOA is not regularly included as a prognostic aerosol species in aerosol microphysics models (including UKCA). Therefore, although the PMOA method is likely the more appropriate representation, we conclude that the SSA method is a reasonable pragmatic approximation for marine-sourced INPs. We note that PMOA emissions are likely more sensitive to future climate and/or inter-annual variability in oceanic productivity than the SSA.

Our results suggest several directions for the improvement of global INP models. There is a clear need for evaluations to be made on longer timescales that can capture diurnal cycles, synoptic-scale events, seasonal cycles, and inter-annual variability. This will provide a robust test of whether models can capture the key transport pathways between source and receptor regions. The evaluations will additionally enable better separation of INPs into particular particle compositions emitted from particular locations, which is needed to build a model that represents sources, transport, and removal of INPs. For this to be possible we require long-term INP measurements that span these timescales. Remote-sensing platforms provide an excellent opportunity to make such measurements, yet complex cloud dynamical and microphysical processes may obscure any link between INP availability and the resulting bulk-scale cloud properties (Korolev et al., 2017; Morrison et al., 2020). This uncertainty associated with cloud microphysical processes restricts the use of remote-sensing observations for quantification of INP concentrations (e.g. Carlsen and David, 2022, and Zhang et al., 2012). An alternative direction is extending current INP measurements to provide greater spatial coverage and longer sampling periods. Several recent campaigns in the Arctic (Sze et al., 2023; Creamean et al., 2022) and those using the Portable Ice Nucleation Experiment (PINE) cloud chamber instruments (Möhler et al., 2021; Brasseur et al., 2022; Adams et al., 2021; Wilbourn et al., 2024; Vogel et al., 2024) are capable of making these measurements. In future work we aim to repeat this evaluation and focus on these longer records, and we encourage continued effort in extending these types of measurements. The analysis of simulated INP concentrations at ambient conditions in Fig. 5 highlights that INPs are active within the lowest 5 km of the atmosphere north of 40° N and south of 40° S. Therefore, we recommend in situ INP measurements from within the lower atmosphere are targeted within these latitudes to gain a representative dataset for evaluation. We also simulate INPs in the tropical mid-troposphere (∼10 km) that may influence tropical high-cloud radiative properties, yet there is a distinct lack of INP measurements at these altitudes. We recommend targeting these altitudes for sampling of air and INP analysis. The in situ sampling should be complemented by building on our existing knowledge of INP activity from surface samples at important source regions, such as Africa, the Middle East, and Australia, and HLD sources that may play a disproportionately important role. Finally, if dust emissions are associated with fertile soils it is likely that marine-sourced INP concentrations in the Southern Hemisphere will become less important, to an extent that it may be reasonable to assume global INP activity can be represented by dust alone. This would simplify the requirements for adding interactive INPs into climate models. This hypothesis would require additional measurements in the SH, especially on longer timescales that can capture the seasonal cycles of dust, SSA, and PMOA and help constrain the primary drivers of INP activity in the Southern Ocean.

The INP distribution we present here provides a valuable tool for the community. Running a fully coupled aerosol microphysics scheme adds considerable computational expense to simulations; hence, cloud-ice microphysics schemes are not commonly coupled to aerosol fields. For simplicity, temperature-dependent parameterizations are often used to calculate INP concentrations (e.g. Meyers et al., 1992, and Cooper, 1986); however, they cannot represent spatial or temporal variability. Others are more complex and are dependent on concentrations of large (over 0.5 µm diameter) non-SSA particles (DeMott et al., 2010) or mineral dust particles (DeMott et al., 2015). Our model builds on these complex parameterizations and provides a 4D climatology of INP concentrations for use in global and regional modelling studies, as well as other relevant applications. Though we do not fully understand the mechanisms involved, the INP distribution using the fertile soil parameterization given by O'Sullivan et al. (2014) is able to best represent the INP measurements. This informs us on the direction of future model developments, with a need to represent dust aerosol emissions from both desert and fertile soils, rather than just a single class.

Appendix A: Details of INP measurements

DeMott et al. (2010)DeMott et al. (2010)DeMott et al. (2010)DeMott et al. (2010)DeMott et al. (2010)Ardon-Dryer and Levin (2014)DeMott et al. (2016)DeMott et al. (2016)DeMott et al. (2016)DeMott et al. (2016)Mason et al. (2016)Mason et al. (2016)Mason et al. (2016)Mason et al. (2016)Mason et al. (2016)Mason et al. (2016)McCluskey et al. (2018b)Price et al. (2018)Si et al. (2018)Irish et al. (2019)Si et al. (2019)Wex et al. (2019)Wex et al. (2019)Wex et al. (2019)Wex et al. (2019)Jiang et al. (2020)Tobo et al. (2020)Sanchez-Marroquin et al. (2020)Sanchez-Marroquin et al. (2021)Tatzelt et al. (2021)Harrison et al. (2022)

Table A1 Details of the INP measurements used for evaluation. The campaign name is included alongside the reference where available. The spatial distribution of measurements can be found in Fig. 6 and individual datasets compared against the INP model in Fig. A1.

Figure A1 Evaluation of the INP model against individual measurement datasets, detailed in Table A1. Dust INPs are represented using the K-feldspar-based parameterization of Harrison et al. (2019). Marine-sourced INPs are represented with the parameterization of McCluskey et al. (2018a) applied to the surface area of SSA. The simulated concentration shown is the monthly median for the corresponding month of the measurement. The solid line corresponds to a bias of zero and the dashed lines to a bias of 1.5 orders of magnitude. Geographical locations, shown in Fig. 6, are depicted by different marker shape/colour combinations.