Quantifying Antarctic‐Wide Ice‐Shelf Surface Melt Volume Using Microwave and Firn Model Data: 1980 to 2021

Antarctic ice‐shelf stability is threatened by surface melt, which has been implicated in several ice‐shelf collapse events over recent decades. Here, we first analyze cumulative days of wet snow/ice status (“melt days”) for melt seasons from 1980 to 2021 over Antarctica's ice shelves using passive and active microwave satellite observations. As these observations do not directly reveal meltwater volumes, we calculate these using the physics‐based multi‐layer snow model SNOWPACK, driven by the global climate‐reanalysis model Modern‐Era Retrospective analysis for Research and Applications Version 2. We find a strong non‐linear relationship between melt days and meltwater production volume. SNOWPACK's calculation of melt days shows agreement with observations of both cumulative days, and spatial and interannual variability. Highest melt rates are found on the Peninsula ice shelves, particularly in the 1992/1993 and 1994/1995 austral summers. Over all ice shelves, SNOWPACK calculates a small, but significant, decreasing trend in both annual melt days and meltwater production volume over the 41 years.

Understanding current and predicting future Antarctic-wide ice-shelf vulnerability to surface meltwater-induced breakup requires a detailed, historical record of meltwater production. Spaceborne microwave sensors and climate models including a detailed firn representation are currently some of the primary tools for analyzing the presence of surface and near-surface melt in Antarctica Johnson et al., 2021;Kittel et al., 2018Kittel et al., , 2021Picard et al., 2022;van Wessem et al., 2018). However, the ability to quantify melt production from these tools is challenging because: (a) microwave sensors can only detect if meltwater is present or not, but they cannot directly quantify meltwater (Picard et al., 2007(Picard et al., , 2022, and (b) climate models generally lack sufficient in-situ snow wetness and density data for validation of their firn components.
Here, we address the limitations outlined above by combining microwave radiometer (i.e., passive) and scatterometer (i.e., active) data with output from the physics-based, detailed, multi-layer snow model SNOWPACK, driven by the global climate reanalysis model Modern-Era Retrospective analysis for Research and Applications Version 2 (MERRA-2) (Gelaro et al., 2017). We first evaluate model performance by comparing annual cumulative melt days (i.e., days when meltwater is present) from SNOWPACK against those derived from the microwave data. Then, we analyze SNOWPACK-calculated meltwater production volume over all ice shelves from 1980 to 2021. We build upon a previous analysis of RACMO2 simulations for the period 1979 to 2010 by Kuipers Munneke et al. (2012), who found an insignificant decrease in surface meltwater production volume, by extending the analyzed time period to 2021, and using a more sophisticated firn model.

Ice Shelves and Ice Shelf Regions
We define ice-shelf boundaries (i.e., calving fronts, grounding lines and margins) using the MEaSUREs v2 data set (Mouginot et al., 2017), which was updated most recently in 2016. To maintain a constant area of the ice shelves through the 41 years of the study, we use these outlines for the entire period, thereby excluding any areas that may have either broken-up or grown in size between 1980 and 2016. For regional melt pattern analysis, we divide the ice shelves into eight regions following Alley et al. (2018) (Figure 1a).

Large-Scale Microwave Radiometer Observations of Melt
Microwave radiometers (19 GHz) record brightness temperatures that depend primarily on the snow temperature and emissivity (Colliander et al., 2022;Zwally, 1977). When liquid water exists in the snow, absorption significantly increases, resulting in increased microwave emissivity, and thus substantially higher brightness temperatures (Liu et al., 2006;Zwally & Fiegles, 1994). Here, we used the near-daily polar stereographic 25 km melt product (version 2) of Torinesi et al. (2003) and Picard and Fily (2006) updated for austral summers from 1980 to 2021. This melt product, widely used before (e.g., Banwell et al., 2021;Magand et al., 2008;Wille et al., 2019), provides binary wet/dry status using observations from five successive sensors from 1980 to 2021 (Scanning Multi-channel Microwave Radiometer and a series of Special Sensor Microwave Imager (SSM/I) sensors). Based on radiative transfer simulations, radiometer brightness temperatures are typically sensitive to the dry/wet status in the uppermost ∼2 m of snow/ice (Picard et al., 2022).
We used data from melt seasons (November 1 to March 31 inclusive) from 1980/1981 to 2020/2021, apart from 1987/1988 due to 41 days of missing data. For each melt season and each 25 km grid cell within the ice-shelf mask, we calculated a daily time series of radiometer-derived wet/dry status, resulting in a calculation of cumulative "melt days". Throughout this paper, we use the phrase "melt days" when referring to days of "wet status" in the microwave data, but note that active melting may not be occurring; meltwater may only be present, and/or may be in the process of refreezing.

Small-Scale Microwave Scatterometer Observations of Melt
We also derived smaller-scale (4.45 km) melt information from C-band (5.225 GHz) radar backscatter images collected by EUMETSAT's Advanced SCATterometer (ASCAT). The 4.45 km enhanced-resolution product was obtained by applying the Scatterometer Image Reconstruction (SIR) algorithm with filtering (Lindsley & Long, 2016) to improve the spatial resolution (Early & Long, 2001). For each day and grid cell, meltwater is assumed to be present when the ASCAT signal is lower than the winter mean signal minus 3 dB (Ashcraft & Long, 2006;Banwell et al., 2021) In locations where snow and firn layers are devoid of liquid meltwater, ASCATS's penetration depth is of the order of meters to tens of meters, hence ASCAT is more likely to be sensitive to water present at depth than microwave radiometers (Picard et al., 2022).
For melt seasons from 2007/2008 to 2019/2020, and for each grid cell within the ice-shelf mask, we calculated the daily time series of scatterometer-derived wet/dry snow status, and hence cumulative melt days.

Firn Model SNOWPACK
To calculate days when meltwater is present in the firn, meltwater production volume, and the presence of surface meltwater versus meltwater at depth, we use the detailed, physics-based, multi-layer model SNOWPACK (Bartelt & Lehning, 2002;Lehning et al., 2002aLehning et al., , 2002b. SNOWPACK has been extensively validated for liquid water flow processes (e.g., Wever et al., 2014Wever et al., , 2015Wever et al., , 2016, and has been used in polar regions including Antarctica

of 11
SNOWPACK was forced with hourly air temperature, relative humidity, wind speed, incoming short-and longwave radiation, and precipitation from MERRA-2 (Gelaro et al., 2017). Gossart et al. (2019) report generally good agreement between MERRA-2 and observed near-surface climatology, including 2 m air temperature and wind speeds, at elevations <500 m in Antarctica, including the ice shelves.
This study analyzes melt in all austral summers from 1980/1981, the first complete summer when MERRA-2 data are available, to 2020/2021. At the ±70° latitude bands, MERRA-2 has a resolution of 24 km × 56 km. To minimize model run time, MERRA-2 points (i.e., grid cell centers) were subsampled based on the similarity in the spatial distribution of climatological conditions following Smith et al. (2020), whereby points were retained in areas of high spatially varying climatic conditions, but removed where climatic conditions were similar to neighboring points. From this sub-selection of MERRA-2 points, we ran SNOWPACK at all points that were located both: (a) on ice shelves >625 km 2 (i.e., at least one 25 km passive microwave grid cell); and (b) <200 m elevation above local sea level according to MERRA-2's model topography. In total, there were 847 MERRA-2 points that fulfilled both criteria.
Using the model setup described by Keenan et al. (2021) and the spin-up procedure from Thompson-Munson et al. (2023), SNOWPACK was first repeatedly run at each of these 847 MERRA-2 points over the period 1980 to 2021 until the spin-ups reached either (a) a >150 m thick column of snow, firn, and ice, or (b) the bottom 2 m of the simulated firn column consisted of solid ice with a total firn column >10 m depth. At 11 grid points, neither of these conditions were reached due to high ablation, restricting the final analysis to 836 grid points, located on 36 ice shelves ( Figure S1 in Supporting Information S1). As the ice-shelf floats on ocean water, we prescribed a fixed lower-boundary temperature of −1.8°C; the freezing point of ocean water. This deviates from Keenan et al. (2021), where the annual average air temperature was prescribed for this lower boundary. While that approach is adequate for dry, grounded ice, it no longer holds for floating ice with surface melt. SNOWPACK calculates meltwater production via a physically based energy-balance model that captures processes such as changes in accumulation, snow-albedo feedback, percolation and latent heat release through the entire firn column (e.g., Wever et al., 2014Wever et al., , 2015Wever et al., , 2016. We expect this to better capture complex melt dynamics compared to, for example, the positive-degree day parameterization using air temperature only, as used in the Community Firn Model (e.g., Medley et al., 2022;Vandecrux et al., 2020). To calculate meltwater volumes over all areas within the ice-shelf mask, we use the method of radial basis functions to interpolate SNOWPACK-simulated melt at each of the 836 MERRA-2 points.

SNOWPACK Validation Against Microwave Data
To evaluate SNOWPACK's performance in calculating realistic spatial and temporal variations in melt, we compared modeled and microwave-derived cumulative melt days for each melt season from 1980/1981 to 2020/2021 for the passive microwave data, and from 2007/2008 to 2019/2020 for the active microwave data, at every one of the 836 MERRA-2 forcing grid points. For each point, we compared the closest radiometer and scatterometer ice-shelf grid cells.
From the SNOWPACK output, we calculated cumulative melt days over different ranges of snow depths from the ice surface (top most 0.1 m, to the complete snowpack depth), and for different thresholds of minimum liquid water content (0-2 mm). We then calculated the root mean square error (RMSE) between the observed and modeled number of melt days, for each of the 836 MERRA-2 points over the 41 melt seasons for the passive data and 13 melt seasons for the active data.
As acquisitions of the radiometers and scatterometers are in the local mornings and evenings (Lindsley & Long, 2016;Picard & Fily, 2006), mean statistics from SNOWPACK were calculated from 6 to 10 a.m. and 6 to 10 p.m. local time (defined by grid point longitude), to more adequately compare to the time-constrained melt observations.

Evaluation of Microwave-Derived and SNOWPACK Calculated Melt Days
We first compare cumulative melt days in each melt season from 1981 to 2020 observed by passive versus active microwave data at each of the 836 MERRA-2 grid point locations (Figure 1b, RMSE = 6.8 days). Our results suggest that for areas of low cumulative melt days (i.e., <40), for example, in the Ronne-Filchner and Ross regions (Figure 1b, Figures S2e and S2f in Supporting Information S1), active microwave observes fewer melt days compared to passive microwave. But for areas with more melt days (i.e., >80), for example, the Peninsula (Figure 1b and Figure S2d in Supporting Information S1) and in particular Wilkins, where microwave-derived melt days are higher than for any other ice shelf ( Figure S3 in Supporting Information S1), the active microwave data reports more melt days than the passive data. These findings make sense, as active microwave is less sensitive to weak surface melt than passive microwave, but active microwave can detect meltwater deeper in the snowpack due to its longer wavelength, and so it often detects meltwater much later into the fall (Picard et al., 2022;Weber Hoen & Zebker, 2020).
Second, we compare cumulative melt days derived from passive and active microwave to those calculated by SNOWPACK over a range of firn depths, and for a range of minimum meltwater detection thresholds (Table S1 in Supporting Information S1) as described in Section 2.5. For both the passive and active microwave data, the SNOWPACK output calculated from the surface to a 0.1 m firn depth and a minimum meltwater threshold of 0 mm produces the lowest mean RMSEs (8.2 and 10.9 days for passive and active, respectively). Using this optimum SNOWPACK output, the bias between the modeled and passive microwave-derived melt days, for all MERRA-2 points and for all melt seasons, is −0.14 melt days (Figure 1c), indicating that SNOWPACK slightly underestimates melt days. Considering the ice-shelf regions separately (Figures 2a-2g), the agreement between modeled and passive microwave melt days is less good. The lowest biases and RMSE values between these data are found in the Ronne-Filchner and Ross Sea regions, despite high scatter there. When just mean melt days in the seven regions are compared, the agreement is again, better (Figure 2h).
Using this same optimum SNOWPACK output, comparing modeled and active microwave-derived melt days for all MERRA-2 points from 2007/2008 to 2021/2022 indicates that, again, SNOWPACK slightly underestimates melt days (Figure 1d; bias = 0.51 melt days). At a regional level ( Figures S4a-S4g in Supporting Information S1), the agreement is, again, less good. However, when the mean melt days in the seven regions are compared, the agreement is again, better ( Figure S4h in Supporting Information S1). Although the precise number of microwave-derived and modeled cumulative melt days differs slightly between the datasets in each region, the patterns of interannual variability are relatively consistent between both microwave sensors and SNOWPACK (Figure 3). Considering all 836 points in all regions (Figure 3a), the RMSE is 68% of the standard deviation, suggesting that the model has some capability in reproducing the observed interannual variability in melt, and also that the climatic conditions are generally well represented in the MERRA-2 reanalysis used to force SNOWPACK.
Within individual regions, we observe local-scale spatial variations in the agreement between passive microwave-derived and modeled melt days ( Figure S5 in Supporting Information S1). For example, SNOW-PACK underestimates melt days in the southwest Peninsula, particularly on the Wilkins and south George VI (bias <−30 days), whereas it overestimates melt days on the Scar Inlet and north Larsen C (bias <15 days). These spatial variations may be partially attributable to the relatively coarse grid of MERRA-2, meaning that SNOWPACK fails to capture localized atmospheric processes (see Section 3.2). Lateral meltwater flow (e.g., Dell et al., 2020) may also contribute to discrepancies.

Modeled Meltwater Volume on Ice Shelves, 1980 to 2021
Over the 41-year time period from 1980/1981 to 2020/2021, modeled mean meltwater volume over all ice shelves is 43.2 Gt yr −1 with high interannual variability (std. dev. = 18.3 Gt yr −1 ) (Figure 4a). Record high and low meltwater volumes were produced in 1997/1998(88.1 Gt) and 1985/1986, respectively. SNOWPACK-calculated mean meltwater volumes over all ice shelves from 1980/1981 to 2020/2021 are lower than the equivalent values from two polar regional climate models, modèle atmosphérique régional (MAR) v3.12  and the regional atmospheric model (RACMO) RACMO2.3p2 (van Wessem et al., 2018(van Wessem et al., , 2023; 74.5 Gt yr −1 and 68.1 Gt yr −1 respectively. However, when we compare SNOWPACK melt (mm w.e.) simulated at the 836 MERRA-2 points to melt calculated in the nearest MAR and RACMO grid cells, for each region and over the 41-year time period, we show that SNOWPACK's melt underestimation is less substantial ( Figure S6 in Supporting Information S1). For MAR versus SNOWPACK, the regression line (R 2 = 0.98) slope is 0.66, and for RACMO versus SNOWPACK, the regression line (R 2 = 0.94) slope is 0.93. In  comparison, for RACMO versus MAR, the regression line (R 2 = 0.94) slope is 1.39, indicating that at these 836 points, modeled melt from the three models is comparable.
As mentioned earlier, differences in the spatial resolution between the higher-resolution models RACMO and MAR and the lower-resolution product MERRA-2 may explain why SNOWPACK melt rates are generally slightly lower. MERRA-2's coarse grid will not be able to capture localized atmospheric processes that often act to increase local melt rates. Such processes include periodic foehn winds over Larsen C (Elvidge et al., 2016), and katabatic winds on Roi Baudouin (Dunmire et al., 2020;Lenaerts et al., 2017) that blow away snow revealing blue ice. For north George VI, the MERRA-2 climate simply cannot be resolved due to high topography on either side of this narrow ice shelf. Therefore we cannot run SNOWPACK here, so no melt is modeled here. As we interpolate point-simulated melt over all ice-shelf areas to calculate volumes ( Figure S7 in Supporting Information S1), melt volumes over north George VI are therefore underestimated compared to the high melt observed here (e.g., Banwell et al., 2021;Holt et al., 2013).
We find a significant (p < 0.05) decreasing trend in SNOWPACK-calculated annual meltwater production volume of −0.55 Gt yr −2 from 1980/1981 to 2020/2021 (Figure 4a). This finding is in line with significant (p < 0.05) decreasing trends in modeled and passive microwave-derived annual melt days (−0.078 and −0.061 d yr −2 , respectively) over the same time period. Additionally, the MERRA-2 mean summer air temperature trend over all ice shelves and over the same time period is negative (−0.016 C yr −2 ), however it is only significant at p < 0.1 (Figure 4a). RACMO2.3 and MAR v3.12 melt volumes over all ice shelves and the same time period (not shown) also show decreasing, albeit not significant, trends of −0.27 and −0.09 Gt yr −2 respectively.
Considering meltwater production volumes in the eight ice-shelf regions (Figure 4b), by far most of Antarctica's meltwater is produced on the Peninsula (mean = 24.6 Gt yr −1 ) where mean melt rates are <435 mm w.e. yr −1 on the Scar Inlet and north Larsen C Ice Shelf ( Figure S7 in Supporting Information S1). After the Peninsula, the next highest meltwater volumes are produced in Wilkes Land (mean 4.5 Gt yr −1 ), the Amery Region (mean 4.4 Gt yr −1 ), and Dronning Maud Land (mean 4.2 Gt yr −1 ) ( Table S2 in Supporting Information S1). Ronne-Filchner, Ross Sea, and Victoria Land all have mean melt volumes <2 Gt yr −1 . Over the 41-year period, all eight regions experience overall decreasing meltwater production rates (Table S2 in Supporting Information S1), though these trends are only significant at p < 0.05 for the Amery, Amundsen Sea, and Dronning Maud Land regions (Wilkes Land is significant at p < 0.1). Notably, the Peninsula's annual meltwater production volume reaches a maximum of 55.6 Gt in 1992/1993 before decreasing at a statistically significant rate (−0.88 Gt yr −2 , p < 0.05) until 2020/2021 (Figure 4b and Figure S7 in Supporting Information S1). We also see a significant decrease in MERRA-2 mean summer air temperature over the same time period (−0.06°C yr −2 , p < 0.05; Figure  S8 in Supporting Information S1).
We find a strong non-linear positive relationship between melt days and meltwater production volume, varying between 0.1 Gt melt day −1 for areas with meltwater production <5 Gt yr −1 (e.g., Ross Sea, Ronne-Filchner) to 1.5 Gt melt day −1 (Figure 4c) for areas with meltwater production >40 Gt yr −1 (e.g., Peninsula in high melt seasons). This finding indicates how microwave data alone, which enables the quantification of "melt days," is not a good indicator of melt volumes. This is because, in warmer summers (Figure 4d), the melting point (0 C) is passed more frequently, the positive temperature-albedo feedback mechanism (Trusel et al., 2015) results in further melting, and larger volumes of meltwater will take longer to refreeze than smaller volumes. Hence we argue for the use of sophisticated energy-balance models such as SNOWPACK.

Conclusions
This study first provided a historical record of cumulative melt days over Antarctica's ice shelves for melt seasons from 1980 to 2021 using passive and active microwave satellite observations. Second, as these microwave observations do not reveal meltwater volumes, these were calculated using the snow model SNOWPACK, driven by MERRA-2. We found a strong non-linear positive relationship between melt days and meltwater production 10.1029/2023GL102744 9 of 11 volume. Evaluation of SNOWPACK's performance in calculating melt days showed agreement with observations in terms of both cumulative days and variability, which gave us confidence in SNOWPACK's ability to also simulate realistic temporal and spatial relative variations in meltwater volumes.
Over the 41-year historical period, SNOWPACK-calculated meltwater volumes were highest on the Peninsula. Over all ice shelves, SNOWPACK showed small, but statistically significant (p < 0.05), decreasing trends in both annual melt days and meltwater production over this time period, in line with a significant decreasing trend passive microwave-derived annual melt days. RACMO2.3 and MAR melt data also showed decreasing, albeit not significant, trends in ice-shelf annual meltwater production over this same time period. As refreezing of meltwater in firn decreases FAC, thereby increasing the potential for surface meltwater ponding and hydrofracture (e.g., Kuipers Munneke et al., 2014), we suggest that the small changes in meltwater production over the past four decades found in this study may have helped to mitigate against ice-shelf instability. Finally, we note that despite the small changes in ice-shelf surface melt to date, projected atmospheric warming means that surface meltwater production on ice shelves is expected to increase non-linearly in the future, and hence ice shelves are predicted to become more vulnerable to future surface meltwater-induced instability Lai et al., 2020). As an example, Gilbert and Kittel (2021) have predicted that 34% of all Antarctic ice shelf areas may be vulnerable to collapse at 4°C of warming above pre-industrial levels. Award #1841607 to the University of Colorado Boulder, and from a University of Colorado Boulder Cooperative Institute for Research in Environmental Studies (CIRES) Innovative Research Proposal (IRP) award. NW was supported by the National Aeronautics and Space Administration (NASA) under Award #80NSSC20K0969 issued by the "Studies with ICESat-2" program. DD was supported by a NASA FINESST Fellowship (Award #80NSSC19K1329). GP received support from the ESA Project ESA/AO/1-9570/18/I-DT-"4D Antarctica." The authors thank Ludovic Brucker, Ted Scambos, Megan Thompson Munson and Rebecca Dell for useful discussions. The authors thank the Editor, Mathieu Morlighem, and two reviewers; Peter Kuipers Munneke and one anonymous reviewer, for their very useful comments.