2.1 In Situ Data Sets and Instrumentation

This study compiled observation data sets of 11 U.S. National Science Foundation (NSF) and three U.S. Department of Energy (DOE) airborne campaigns, including NSF START08 (Pan et al., 2010), HIPPO (Wofsy, 2011), PREDICT (Montgomery et al., 2012), TORERO (Volkamer et al., 2015), DC3 (Barth et al., 2015), CONTRAST (Pan et al., 2017), WINTER (Lee et al., 2018), CSET (Albrecht et al., 2019), ORCAS (Stephens et al., 2018), SOCRATES (McFarquhar et al., 2021), and OTREC (Fuchs-Stone et al., 2020) campaigns, as well as DOE M-PACE (Verlinde et al., 2007), ISDAC (McFarquhar et al., 2011), and ACME-V (Maahn et al., 2017) campaigns. Table S1 in Supporting Information S1 provides detailed information about these campaigns, including name, number of research flights, time, location, and flight hours at all temperatures as well as at mixed-phase range. A total of 576 hr were flown at −35°C–0°C, with 463 and 113 hr of clear-sky and in-cloud conditions, respectively. Collectively, the 14 airborne campaigns conducted between April 2008 and September 2019 covered a spatial range from 75°S to 87°N and from 38°W to 128°E (Figure 1a). A similar but smaller in situ observation data set was previously used in Wang et al. (2024) as described in their Table 1 to validate three satellite-based cloud phase products.

For most campaigns (except M-PACE and ISDAC), we applied a phase identification method that involves (a) the identification of ice or supercooled liquid for each second of a cloud probe and (b) IWC and LWC calculations (D’Alessandro et al., 2019; Maciel et al., 2024; Yang et al., 2021). Several key parameters were used, including the mass–number concentration (M–N) relationship, maximum particle diameters, the standard deviations of particle size distributions, and temperature. For most campaigns, these parameters were derived from the Cloud Droplet Probe (CDP) (2–50 μm) and the Fast-Two Dimensional Cloud probe (Fast-2DC) (62.5–3,200 μm). The 2-Dimensional Stereo (2DS) cloud probe (40–5,000 μm) was used instead of the Fast-2DC probe in SOCRATES and ACME-V. The 1-Hz CDP observations were categorized as large aerosols (N ≤ 10^−1.5 cm⁻³ or M ≤ 10^−3.4 g m⁻³), ice crystals (both 10^−1.5 < N < 10^−0.5 cm⁻³ and M > 10^−3.4 g m⁻³), or liquid droplets (both N ≥ 10^−0.5 cm⁻³ and M > 10^−3.4 g m⁻³). A more complex decision tree was applied for Fast-2DC and 2DS following D’Alessandro et al. (2019).

To calculate LWC or IWC, spherical shape was assumed for supercooled liquid droplets, while the mass-dimension (M-D) relationships from Brown and Francis (1995) were used for small and large ice particles, separated by maximum dimensions ≤ and >75 μm, respectively. M-PACE (McFarquhar et al., 2007) and ISDAC (Jackson et al., 2012) derived IWC and LWC from two different suites of instruments. M-PACE used Forward Scattering Spectrometer Probe (FSSP), 1-Dimensional Cloud probe (1DC), 2DC, and High Volume Precipitation Spectrometer (HVPS), while ISDAC used FSSP, CDP, 2DC, 2DS, Cloud Imaging Probe 2 (CIP2), and 2-D Precipitation (2DP) probe.

A consistent definition of three cloud phases was applied to all observations and simulations. That is, if supercooled liquid fraction (SLF) < 0.1, between 0.1–0.9, or >0.9, then this sample is defined as ice, mixed, or liquid phase, respectively. Here SLF = LWC/(LWC + IWC), which represents the mass fraction of supercooled liquid water.

2.2 Model Simulations

The DOE Energy Exascale Earth System Models version 1 (E3SMv1), specifically its atmospheric component, the E3SM Atmosphere Model version 1 (EAMv1), was used in this work (Rasch et al., 2019; Xie et al., 2018). EAMv1 employs a Spectral Element (SE) dynamical core at ∼1-degree resolution with 72 vertical levels (Golaz et al., 2019). The SE dynamical core uses unstructured grids and has advantages of near-perfect scalability and GPU (Graphics Processing Unit) acceleration (Adbi et al., 2017; Dennis et al., 2012). Cloud microphysics in EAMv1 are treated by version 2 of Morrison and Gettelman (MG2) with representations of cloud particle formation, growth, and precipitation processes (Gettelman & Morrison, 2015). Cloud macrophysics, shallow convection, and boundary layer turbulence are simulated using the Cloud Layers Unified By Binormals (CLUBB) (Golaz et al., 2002; Larson & Golaz, 2005), deep convection processes are parameterized based on G. J. Zhang and McFarlane (1995), and the aerosol module used the Modal Aerosol Module (MAM4) (Liu et al., 2016).

For each flight campaign, one model simulation was set up to cover the same time period with ∼6 months of spin-up time prior to the start date of the campaign. The model output was saved along the flight tracks at a 10-min frequency. Simulations were nudged toward the Modern-Era Retrospective Analysis for Research and Applications, Version 2 (MERRA-2) reanalysis data set (Gelaro et al., 2017) for temperature and horizontal wind fields, consistent with previous model evaluation studies (e.g., Desai et al., 2023, 2025; Patnaude et al., 2024; Yang et al., 2021; Yip et al., 2021). Several sensitivity tests were conducted to examine the impacts of vertical sampling and horizontal spatial averaging, similar to the method used in Yang et al. (2021) in their Section 2.3, and consistent results were found among these tests.

3.1 Cloud Phase Frequency Distributions and Phase Partition at a Near Global Scale

One advantage of airborne observations is the ability to provide synchronized in situ measurements of cloud hydrometeors and environmental conditions (such as temperature). Hence, the mass partition between supercooled liquid and ice (indicated by SLF) is examined as a function of temperature in two hemispheres (Figures 1b and 1c), with the number of samples shown in Figure S1 of Supporting Information S1. As for seasonal distributions, most samples in the SH occurred during austral summer, while the NH has more similar sample sizes in all four seasons (Table S2 and Figure S2 in Supporting Information S1).

Based on high-resolution 1-s observations, a main hemispheric difference is the higher SLF in the SH compared with the NH, with SLF in the SH being 0.2–0.4 higher at −25 to 0°C. A moving average at 100-km resolution that included both clear-sky and in-cloud segments was used to compare with model grid-mean values, consistent with the method used in D’Alessandro et al. (2019) and Yang et al. (2021). After spatial averaging, SLF shows larger increases in the NH (by 0.1–0.2) than the SH (by 0.05–0.1), because the 1-s observations are dominated by pure ice segments especially in the NH, while 100-km observations have more mixtures of ice and liquid.

Compared with the 100-km observations, the EAMv1 simulations show similar SLF values in the SH but show much higher SLF in the NH by 0.1–0.2. Consequently, the EAMv1 simulations do not show a significant contrast of SLF between the two hemispheres, while 100-km observations show higher SLF in the SH than NH. This result indicates that future improvement of the model parameterization should convert more liquid-containing clouds to the ice-containing clouds at temperatures between −35°C and 0°C in the NH.

3.2 Latitudinal-Temperature Distributions of Three Cloud Phases

The latitudinal-temperature distributions of three cloud phases are shown in Figure 2. The number of samples and in-cloud frequencies are shown in Figures S3 and S4 of Supporting Information S1, respectively. To reduce the noise in frequency distributions due to fluctuations of sampling sizes, we applied an averaging process to every 3 × 3 bins (i.e., a center bin and its 8 surrounding bins). Key features are consistently seen without the 3 × 3 grid smoothing (Figure S5 in Supporting Information S1) or using a larger model sample size (Figure S6 in Supporting Information S1). The cloud phase frequency is calculated as the number of a phase divided by the total number of in-cloud samples in each bin. In addition, in-cloud frequency is calculated as the number of in-cloud samples divided by the total number of all-sky samples.

A hemispheric contrast is consistently seen in this latitudinal view based on in situ observations at both ∼0.2 and 100-km resolutions, showing higher liquid phase frequency in the SH compared with the NH (Figures 2a and 2b). The EAMv1 simulations show similar latitudinal gradient of in-cloud frequencies compared with 100-km observations, with the highest in-cloud frequencies located at the polar regions (Figure S4 in Supporting Information S1). However, the lack of hemispheric differences in frequency distributions of liquid and ice phases is seen in the simulations, consistent with the lack of hemispheric asymmetry in SLF (Figure 1).

The hemispheric differences for observations and simulations are more evidently represented in Figure 3, including the differences between each pair of latitudinal-temperature bins symmetrically distributed between the NH and SH (columns 1–4) and the phase occurrence frequencies in each 5°C temperature bin (column 4). The number of samples are shown in Figure S7 of Supporting Information S1. Both 1-s and 100-km observations show higher liquid phase frequency in the SH by 0.2–0.4 than the NH at −25 to 0°C (Figure 3d). Such observed hemispheric differences are even larger at higher temperatures than lower temperatures. In addition, both observations show significantly higher ice phase frequencies in the extratropical regions in the NH (>30°N) compared with the extratropics in the SH (<30°S) (Figures 3i and 3j).

The model simulations show lack of hemispheric differences in both liquid and ice phase frequencies (Figures 3d and 3l). These lack of hemispheric differences in the simulations can be attributed to the different magnitudes of model biases between the two hemispheres for liquid and ice frequencies. That is, simulations in the SH show similar liquid and ice phase frequencies compared with 100-km observations (Figures 3h and 3l). However, in the NH, the model underestimates ice phase frequency and overestimates liquid phase frequencies from −30 to 0°C.

Smaller model biases are seen for mixed phase frequencies when compared with the 100-km observations, especially for the SH (Figure 3h). However, the model still overestimates mixed phase frequency at −30 to −15°C and underestimates it at −15 to 0°C in the NH. In addition, the model misses the trend of increasing mixed phase frequencies at higher temperatures as shown by the observations. The results from Figures 2 and 3 indicate that the model biases in terms of phase occurrence frequency are more severe in the NH, consistent with the SLF biases in Figure 1.

3.3 Hemispheric Comparisons of Cloud Microphysical Properties in Respective Cloud Phases

Cloud microphysical properties, that is, LWC and IWC, are examined for individual cloud phase or all phases at various temperatures (Figures 4a–4h) alongside their hemispheric differences (i.e., NH minus SH) in Figures 4i–4p. The number of the samples is shown in Figure S8 of Supporting Information S1. Due to the definitions of three cloud phases, each phase may contain a certain amount of IWC, LWC, or both.

Focusing on the dominate type of hydrometeors in the liquid and ice phases, both observations show higher LWC of the liquid phase in the SH than NH (Figure 4a) and relatively similar IWC of ice phase between the two hemispheres (Figure 4g). Spatial averaging leads to more pronounced decreases of LWC in liquid phase (Figure 4a) and smaller decreases of IWC in ice phase (Figure 4g) because most of the liquid segments are shorter than ice segments (not shown). More significant decreases of LWC are seen in the NH than SH after averaging, consistent with the more extensive coverage of supercooled liquid clouds over the SH as reported by previous studies (e.g., Barone et al., 2024; Desai et al., 2023, 2025; Hu et al., 2010).

Compared with 100-km observations, EAMv1 significantly overestimates the LWC in the liquid phase (Figure 4a) and mixed phase (Figure 4b), as well as overestimating the total LWC of all phases (Figure 4d) by 1–2 orders of magnitude. In fact, the simulated LWC is closer to the 1-s observations than the 100-km observations. The overestimations of simulated LWC in both hemispheres lead to small hemispheric differences of LWC at −25 to −10°C, similar to the 100-km observations at that temperature range (Figure 4l). At temperatures above −10°C, both observations show higher total LWC in the SH than NH, but the model shows the opposite hemispheric difference.

Differing from the consistent overestimation of LWC in all phases, the model biases in IWC vary with phase. The simulated total IWC of all phases is more similar to the 100-km observations except for the large negative biases in the NH around −20 to −5°C (Figure 4h). However, the simulations overestimate IWC for mixed phase in two hemispheres (Figure 4f) and underestimate IWC for ice phase in the SH (Figure 4g). Despite these IWC biases for individual phase, the hemispheric differences in the total IWC (Figure 4p) show similar results between EAMv1 and 100-km observations at −30 to −20°C with higher total IWC in the NH. However, EAMv1 misses the higher total IWC in the NH at −20 to 0°C due to its underestimation of the total IWC in the NH at that temperature range (Figure 4h).