Drivers of ASCAT C band backscatter variability in the dry snow zone of Antarctica

Abstract

C band backscatter parameters contain information about the upper snowpack/firn in the dry snow zone. The wide incidence angle diversity of the Advanced Scatterometer (ASCAT) gives unprecedented characterisation of backscatter anisotropy, revealing the backscatter response to climatic forcing. The A (isotropic component) and M₂ (bi-sinusoidal azimuth anisotropy) parameters are investigated here, in conjunction with data from atmospheric and snowpack models, to identify the backscatter response to surface forcing parameters (wind speed and persistence, precipitation, surface temperature, density and grain size). The long-term mean A parameter is successfully recreated with a regression using these drivers, indicating strong links between the A parameter and precipitation on long timescales. While the ASCAT time series is too short to determine which factors drive observed trends, factors influencing the seasonal and short timescale variability are revealed. On these timescales, A strongly responds to the propagation of surface temperature cycles/anomalies downward through the firn, via direct modulation of the dielectric constant. The influence of precipitation on A is small at shorter timescales. The M₂ parameter is controlled by wind speed and persistence, through modification of monodirectionally-aligned surface roughness. This variability indicates that throughout much of coastal Antarctica, a microwave ‘snapshot’ is generally not representative of longer-term conditions.

1. INTRODUCTION

Satellite-based microwave observations at C band (4 to 8 GHz) have been used for ice sheet near-surface parameter retrieval for several decades (e.g. Spiridonov and Ozerkina, 1984; Bingham and Drinkwater, 2000; Drinkwater and others, 2001; Winebrenner and others, 2001; Arthern and others, 2006; Fraser and others, 2014). Observations indicate a dry snow penetration depth (e⁻¹) of ~20 m for C band radar (Rignot, 2002), in agreement with theory (e.g. Bingham and Drinkwater, 2000), and 99% of the backscatter is contributed by the top ~90 m (Partington and Flach, 2003). C band active microwave returns are sensitive to snow grain size, temperature, net annual surface mass balance (SMB) and other parameters. However, with the exception of very few studies (e.g. Bingham and Drinkwater (2000), which covered only north of 72°S because of the instrument orbit), previous parameter retrieval in polar ice sheets has always been conducted using either long-term means of microwave parameters (e.g. Drinkwater and others, 2001) or temporal snapshots (e.g. Dierking and others, 2012), without consideration of short timescale (i.e. monthly or shorter) variability. A quantification of such variability on a variety of timescales, and an exploration of potential drivers, are the main aims of the present study.

The first Advanced Scatterometer (ASCAT; Figa-Saldaña and others, 2002) was launched in 2007 on board the European Organisation for the Exploitation of Meteorological Satellites (EUMETSAT) MetOp-A satellite. ASCAT is a C band fan-beam scatterometer with a beam configuration that provides excellent simultaneous diversity in both azimuth and incidence angles, particularly equator-ward of ~78°N/S (Fraser and others, 2014). ASCAT also has almost complete polar coverage, with accurate parameter retrieval to ~88°S. A second ASCAT instrument was launched in 2012 onboard MetOp-B, with a third to follow in 2018 (Klaes and others, 2013). The work presented here uses data only from ASCAT on MetOp-A.

Multi-angle backscatter observations of the Antarctic upper snowpack generally show a strong anisotropy, both in azimuthal and incidence angle. Retrieval of properties of the upper snowpack thus relies heavily on characterisation and/or removal of this anisotropy. Incidence angle backscatter anisotropy is well modelled as a linear function of incidence angle in most regions of the Antarctic ice sheet (AIS), although several regions (e.g. ice shelves, high on the East Antarctic plateau) show a weakly cubic response (Fraser and others, 2014). Azimuthal anisotropy is primarily generated by the interaction between incident down-welling microwave radiation and monodirectionally-aligned Rayleigh scale roughness of both the surface and internal layer boundaries. When the incident microwave radiation is aligned perpendicular to the direction of surface relief (i.e. more roughness is encountered by incident microwaves), backscatter is enhanced. These microrelief features are formed by and associated with erosional processes responsible for the formation of sastrugi. In many parts of the AIS, particularly in non-coastal sloping regions where katabatic winds predominate, wind direction is very persistent (Parish and Bromwich, 1987), resulting in strongly mono-directional surface relief (Long and Drinkwater, 2000). The presence of such microrelief results in a strong (>5 dB in some regions) bi-sinusoidal modulation of backscatter as a function of azimuth angle (Fraser and others, 2014).

Given this non-negligible anisotropy, it is natural for scatterometer data to be presented by parameterising surface backscatter measurements as a function of incidence angle and azimuth angle: σ ₀ = A + f(θ) + g(ϕ), where σ ₀ is the measured backscatter, θ is the observation incidence angle and ϕ the azimuth angle. Fraser and others (2014) give a detailed comparison of various ASCAT parametrisations for incidence angle and azimuth angle, concluding that a linear incidence angle parametrisation combined with a Fourier series parametrisation for azimuth angle (with cos(1ϕ), cos(2ϕ) and cos (4ϕ) terms) accounted for much of the observed anisotropy of the AIS without over-fitting. This parametrisation is represented by

(1) $$\eqalign{{\sigma _0} = & A + B(\theta - 40^\circ ) + {M_1}\cos (1(\phi - {\phi _1})) \cr & + {M_2}\cos (2(\phi - {\phi _2})) + {M_4}\cos (4(\phi - {\phi _4}))}$$

where A is the isotropic component of backscatter, B is the (always negative) change in backscatter with increasing incidence angle (here normalised to 40°) and M ₁, M ₂, M ₄, ϕ ₁, ϕ ₂ and ϕ ₄ are the Fourier series coefficients describing the azimuthal modulation. M ₂, the magnitude of the bi-sinusoidal term, is generally dominant over M ₁ and M ₄.

The microwave backscatter characteristics of deep snow and firn undergo irreversible changes after incorporation of liquid water into the solid matrix, even after liquid water has refrozen. Scatterometers are particularly sensitive to such changes (Trusel and others, 2012). For this reason, we analyse only those parts of the AIS that experience very low annual melt (Trusel and others, 2012), here thresholded at a value of 0.2 days of melt a⁻¹ (hereafter referred to as the dry snow zone).

The microwave signature of deep, dry snow/firn (i.e. depth > microwave penetration depth, or ~90 m in this case) is complex, and sensitive to many parameters. The reader is directed to Yurchak (2009) (dealing mainly with the volume component), Ulaby and others (1996) and Partington and Flach (2003) for more complete descriptions, but a short, qualitative summary is provided here, based largely upon these references. Surface backscatter for a rough surface depends on dielectric properties, Rayleigh scale roughness and roughness anisotropy (giving rise to the aforementioned M ₂ azimuthal modulation). Volume backscatter at C band is largely controlled by grain size within the penetration depth (low backscatter from small grains; Figure 16 of Yurchak (2009)), temperature (warmer snow gives lower backscatter) and layer thickness. In the AIS interior, grains are deposited at an initial radius on the order of 0.3 mm (Linow and others, 2012). Grain growth in the upper snowpack is driven by snowpack temperature and the vertical temperature gradient, following an Arrhenius-type growth model (e.g. Alley and others, 1982; Flanner and Zender, 2006). This is balanced by the input of fresh, small grains in accumulating regions. Upper snowpack density variations are also of similar phase, peaking in late summer as sintering and metamorphosis processes are enhanced during the summertime months (Ligtenberg and others, 2011).

Bingham and Drinkwater (2000) developed a microwave backscatter model for multi-layered snow/firn. The model takes initial grain size, temperature mean and cycle amplitude, layer thickness (related to net annual SMB), surface density and grain growth rate as inputs and generates a snow stratigraphy. Volume backscatter is modelled as Rayleigh backscatter from spherical scatterers. There is no rough surface backscatter component in the model. This model was applied to time series of Antarctic (Bingham and Drinkwater, 2000) and long-term mean Greenland (Drinkwater and others, 2001) scatterometer observations, in order to perform parameter retrieval. The Antarctic study attributed seasonal A parameter cycles mainly to surface temperature cycles, and interannual trends were attributed to changes in SMB. The Greenland study by Drinkwater and others (2001) empirically linked A and B to long-term mean net SMB, and used the backscatter model to verify these empirical relations. Using the same model, the sensitivity of modelled backscatter to a wide range of these input parameters is covered in detail by Partington and Flach (2003).

Because of the sensitivity of C band backscatter to grain size and surface roughness, a different signature is expected from regions experiencing consistent, widespread (on a sub-pixel level) positive SMB, compared with other SMB regimes. Megadunes are transverse dunes of small (several metre) amplitude but large (~2 km) wavelength (Fahnestock and others, 2000). On the windward side, SMB is positive, and the surface is typically characterised by fine grains and sastrugi. On the more spatially extensive leeward side, SMB is low (typically <20 kg m⁻² a⁻¹) or negative due to wind-driven snow redistribution/removal or ablation due to sublimation (Scambos and others, 2012; Das and others, 2013). Megadune regions exist in the interior of the East Antarctic ice sheet (EAIS), are widespread and have a strong C band emission and backscatter signature (Fahnestock and others, 2000), so must be considered in this work. The low-SMB leeward side is typically referred to as a glaze region, due to the strongly forward-scattering surface reflection. In glaze regions, SMB is so low that the layer burial paradigm dominant in non-glaze regions breaks down, leading to repeated summertime grain growth in the upper few metres. Glazed surfaces/sub-surfaces are thus characterised by very large grain radii (on the order of 1 cm), a shiny appearance (forward-scattering enhancement) caused by scour and sintering processes, and summertime surface cracks caused by thermal expansion (Scambos and others, 2012). Furthermore, glaze surfaces are generally free from significant sastrugi and centimetre-scale roughness features, since there is little new snow to be eroded by winds. Glaze-prevalent regions occupy ~11.2% of the EAIS above an elevation of 1500 m (Scambos and others, 2012), and are detectable in coarse-resolution scatterometer imagery (Fraser and others, 2014). Glaze regions are an example of where a layer-type backscatter model is not expected to perform well. Not only is the layer assumption questionable in a glazed subsurface where SMB is near-zero, but the snow–air interface itself has very different characteristics to that of layered snowpack.

Aside from characterising backscatter variability, another aim of this work is to elucidate the links between atmospheric driving factors, and the response of C band backscatter across the dry snow zone. This is achieved in a variety of ways. Firstly, the long-term mean values of A and M ₂ are recreated using multiple linear regression (MLR), taking input from various atmospheric and surface datasets. The variability of backscatter is then analysed on three timescales: Interannual, to determine if long-term trends are present in the ASCAT data; intraannual, to quantify and investigate the causes of variability throughout the annual backscatter cycle; and 5 days timescale, to assess, quantify and discuss the short timescale variability.

This study is conducted under the assumption that the observed ASCAT parameter variability is driven by changes in the physical properties of the snow surface and upper snowpack, including changes in density, layer roughness, grain size and temperature. It is assumed that changes to the dielectric properties of snow and ice from chemical species variability is negligible, as has been implicitly assumed in prior studies. In turn, these snowpack and surface physical changes are assumed to be driven solely by atmospheric processes: while intraannual changes to the radiation field (particularly shortwave radiation) are profound in Antarctica, these are difficult to accurately quantify because of the difficulties associated with accurately modelling the cloud field in Antarctica (Inoue and others, 2015).

This study makes use of lagged correlation analyses to study the relationship between atmospheric drivers and ASCAT parameter response. In the absence of a radiative transfer model, it is difficult to conclusively attribute the observed variability to changes in the snowpack/surface, and then to relate these changes to atmospheric drivers. Here, knowledge of physical processes is used to guide interpretation, and to assign probable causality. Furthermore, the sign and magnitude of the lag required to maximise the correlation also contains information on the processes responsible for the correlation (e.g. any correlation where the surface dataset leads the backscatter response clearly violates causality). As such, this work should be interpreted as a study of likely drivers of backscatter parameters, rather than a proven causal link. While the Partington and Flach (2003) model may have been of some utility for this purpose, it does not include a rough surface scattering component, so is of limited use (particularly for studies of the M ₂ parameter, which is controlled primarily by surface roughness). Furthermore, since it was developed primarily for use with the Greenland ice sheet, this model assumes a layered subsurface. This assumption likely breaks down in the spatially-extensive glaze-prevalent regions on the East Antarctic plateau. Development of a replacement model (or modification of previous models) is outside the scope of this paper.

2. DATASETS

ASCAT Level 1B Hierarchical Data Format version 5 (HDF5) data (code ASCA_SZR_1B) were acquired from the EUMETSAT Data Centre (http://www.eumetsat.int). These swath-level data were then reprojected onto a 12.5 km resolution polar stereographic projection defined by the National Snow and Ice Data Centre (NSIDC, grid version 3, details available at https://nsidc.org/data/polar_stereo/ps_grids.html). Following reprojection, parametrisation was performed using a linear incidence angle model and a Fourier series with cos(1ϕ), cos(2ϕ) and cos(4ϕ) terms for azimuthal modulation. This parametrisation was shown to effectively describe the observed anisotropy on the AIS, while not over-fitting (Fraser and others, 2014).

From the eight parameters available, A and M ₂ were chosen for analysis because of their relative ease of interpretation, high signal-to-noise ratio and expected relation to physical properties, as detailed in Section 1. Data prior to July 2007 were discarded due to suspected poorer calibration.

Surface temperature, 10 m wind velocity and precipitation data were obtained from the Antarctic Mesoscale Prediction System (AMPS) model, from July 2007 until January 2012 (Powers and others, 2012). A short-term comparison of various operational atmospheric model and reanalysis products showed that AMPS gives realistic precipitation values at Mill Island, located on the East Antarctic coast at ~100.67°E; 66.5°S, with a summit elevation of ~500 m (Fig. 1). Although net surface accumulation is a more glaciologically-relevant parameter than precipitation, the latter was chosen as a co-variate due to difficulties associated with accurate modelling of processes such as blowing snow redistribution and sublimation after re-suspension (i.e. precipitation is implicitly assumed to be a first-order approximation for net surface accumulation throughout the continent). Grid spacing was 20 km prior to October 2008, and 15 km thereafter. An advantage of this high spatial resolution is more accurate representation of the surface topography, leading to better characterisation of katabatic winds (thought to be an important driver of the M ₂ parameter (Long and Drinkwater, 2000)).

Fig. 1. Plot of precipitation from various atmospheric model products vs Mill Island ice core-derived accumulation (solid lines). Since the AMPS data began in March 2007, the 2007 AMPS data point has been increased by a factor of 1.2 for comparison against the other datasets, effectively assuming equal precipitation in all months. For the RACMO2.1 and AMPS datasets, the closest point to the north-east of the ice core location was chosen, because of the orographically-driven nature of snowfall around the East Antarctic coast. ‘ERA-Interim’ is the closest grid point to the Mill Island summit for the 1.5° grid spacing European Centre for Medium-Range Weather Forecasts interim reanalysis dataset (Dee and others, 2011). ‘RACMO2.1’ is the precipitation field from the RACMO2.1 dataset (27 km grid spacing). Since the orography of the Antarctic coast is an important factor determining precipitation, it is not unexpected that the model with the highest spatial resolution (AMPS, 20 km spacing prior to October 2008; 15 km thereafter) produces the best precipitation record in regions of steep orography such as Mill Island.

While the azimuth of the Antarctic katabatic wind has been remotely retrieved on a large scale using the ϕ ₂ parameter from scatterometer data (Remy and others, 1992; Long and Drinkwater, 2000), explicit comparisons between wind speed/persistence and M ₂ value have not been conducted on large scales. However, it stands to reason that a high M ₂ value is only achieved if the wind blows persistently in a single direction. And while the wind field history must play a role in M ₂ enhancement, a more recent wind field likely contributes more strongly to M ₂ than an older one. The timescale of required wind persistence has yet to be studied. Thus, in order to address this, and to find a stronger link between wind and M ₂, we define P, a wind persistence metric:

(2) $$P = \displaystyle{{\sqrt {{{(u{^\ast}w)}^2} + {{(v{^\ast}w)}^2}}} \over {s{^\ast}w}}$$

where u is the mean zonal wind component, v is the mean meridional wind component and s is the mean wind speed. The exponential weighing factor, w, is defined as w = e^{−t*t _s} . Here, t is time, and the timescale t _s, (d), is a free parameter. The t _s value was selected by maximising the non-glaze-prevalent region (i.e. where the correlation between M ₂ and persistence takes a positive value) mean correlation between the persistence metric and the M ₂ variability. A maximum correlation value, r = 0.24 (mean value for the dry snow zone) is achieved at a value of t _s = 1/0.06, ≃ 16.7 days. Although a value of r = 0.24 is a low correlation value, it shows considerable structure, and locally much higher values (up to r = 0.7) occur on coastal slopes and in glacier valleys, where katabatic winds dominate the wind field. This correlation coefficient is higher than the simple 5 days history wind persistence vs M ₂ correlation value (r = 0.18), and also higher than the simple wind speed vs M ₂ correlation value (also r = 0.18). The temporal window for creating P was truncated at 40 days, where w diminishes to a value of 0.1. According to this model, ~50% of contributions to M ₂ enhancement arise from the most recent 10 days wind activity, providing a characteristic timescale for maximum M ₂ enhancement from sastrugi/longitudinal dune/surface microrelief formation.

The uppermost 200 model layers of snow/firn density data (typically extending to a depth of >2 m) were provided by a firn densification model (Ligtenberg and others, 2011) driven by the regional atmospheric climate model RACMO2.1 (Lenaerts and others, 2012). Grid spacing was 27 km throughout the time series. The data were provided at a monthly time step. Since most atmospheric-driven variability is expected to occur in the part of the snow/firn most exposed to atmospheric drivers, and since radar returns exponentially decay with depth, these data were first converted from layer density into an upper 1 m density dataset. Upper 1 m mean snow grain radius from a snow grain size model (Kuipers Munneke and others, 2011), also forced by RACMO2.1 (Lenaerts and others, 2012) was reprojected to match the NSIDC grid. Together with the surface temperature, wind speed, wind persistence and precipitation, these six datasets are hereafter collectively referred to as the surface datasets. These datasets are considered as possible drivers of the scatterometer datasets (the A and M ₂ parameters) in the following sections. The AMPS fields (surface temperature, wind speed, wind persistence and precipitation datasets) have the potential to indirectly drive backscatter changes via changes in the physical properties of the snowpack/firn. On the other hand, the RACMO-derived upper 1 m mean density and grain size datasets are considered as direct potential drivers, since these fields directly represent snow/firn changes. All datasets cover the time period from July 2007 to January 2012.

A scatterometer-derived map of annual mean melt duration (Trusel and others, 2012) was used, together with a threshold of <0.2 days mean annual melt, to define the dry snow zone, since ASCAT A and M ₂ returns were found to be generally insensitive to melt below this level. A map showing mean annual melt duration and the dry snow zone is given in Figure 2. This figure also indicates the names of regions referred to in the text.

Fig. 2. Map of the AIS, showing mean annual melt duration (rainbow colour scale) from the Trusel and others (2012) dataset, sub-sampled to a 50 km px⁻¹ resolution; and the dry snow zone (grey region, defined as mean annual melt duration <0.2 days). The coastline is rasterised from the MODIS Mosaic of Antarctic (MOA) product (Scambos and others, 2007). Regions mentioned in the text are shown as darker grey ellipses.

Slope aspect is known to strongly influence the net mean SMB, particularly in coastal locations (e.g. Law Dome, East Antarctica; Goodwin and others, 1990). This temporally static field is used here only as a covariate in the MLR analysis of long-term means. This field was derived from the Bamber and others (2009) digital elevation map.

3. METHOD AND RESULTS

3.1. Recreation of long-term mean values using regression

Assuming the atmosphere determines the physical state of the snow surface/upper snowpack, it should be possible to recreate the spatial distribution of ASCAT parameters using spatial time series of atmospheric parameters. One simple tool for achieving such a recreation is MLR. The MLR output for each grid cell is described in the equation (for the A parameter):

(3) $${A_{recreated}} = \sum\limits_{i = 1}^n {k_i}{\,p_i}$$

where n is the number of covariate datasets, p _i are the covariates (e.g. wind speed, temperature, wind persistence), and k _i are coefficients calculated during the MLR process. Here, covariates vary spatially across the dry snow zone, and values of k _i are constant.

MLR works well when the following conditions are met: at least one covariate has a strong correlation with the dependent variable (A or M ₂ here) and all relationships between the covariates and the dependent variable are linear. In the case of nonlinear relationships, often a transformation of a covariate can increase the linearity of the relationship, sometimes at a cost of reduced sensitivity to the dependent variables. In this work, interactions between covariates was examined by including products of two covariates. This added complexity but only reduced the root mean square (RMS) residual marginally in each case, so the results are not given here.

MLR was able to effectively recreate the A parameter when using mean precipitation, surface temperature, wind persistence, slope aspect and upper 1 m mean grain size and density as covariates, as shown in Figures 3a–c. The A parameter response to precipitation displayed strong nonlinearity, so a transformation precip^′ = (precip + 0.6)⁻¹ was employed to linearise the relationship. The addition of 0.6 mm w.e. was chosen by maximising the correlation coefficient between precip^′ and A. Likewise, the relationship between A and slope aspect was linearised using the transformation aspect^′ = abs(aspect + 11°), with the 11° offset also chosen by maximising the correlation coefficient with A. The RMS residual was 1.64 dB in this case (distribution of residual given in Fig. 3c). Given the ~20 dB range of A, this represents an RMS error of ~8.2%.

Fig. 3. MLR of covariates slope aspect (d, g), precipitation (e, h), temperature (f, i), wind persistence (j, m), upper 1 m mean grain size (k, m) and upper 1 m mean density (l, o) to recreate A. (a) The long-term mean A parameter. (b) The MLR recreation of A using the aforementioned covariates. (c) The residual (RMS residual = 1.64 dB or ~8.2% of the A range).

The success of MLR in recreating A is largely a consequence of the strong correlation between A and transformed precipitation (Figs 3e and h, r = 0.78). Other covariates act more subtly to reduce the residual. This is regarded as very strong evidence that precipitation is the dominant driver of A on long timescales. Such a hypothesis is supported by the Bingham layered snowpack model (Bingham and Drinkwater, 2000), as demonstrated in Greenland by Drinkwater and others (2001). The favoured physical interpretation is that higher precipitation gives a lower backscatter by increasing the one-way loss factor L(θ ^′) = e ^{k _e d.secθ ′} from a thicker snow layer, where k _e is the reciprocal of the penetration depth, d is the snow layer thickness and θ ^′ is the refracted incidence angle in the snow (Drinkwater, 1989).

Temperature also shows a reasonable negative correlation with A (r = −0.55). The negative correlation implies that this relationship is not due to snow grain growth, since higher temperatures should enhance grain growth leading to higher backscatter. With this mechanism excluded, alternative possible physical interpretations include direct modification of the dielectric constant with temperature (as shown by Partington and Flach, 2003), or a Clausius-Clapeyron-type relation whereby less precipitation falls in colder locations. There is a strong relationship between A and temperature at shorter time-scales (as shown in the sections below), so the dielectric constant modification hypothesis is preferred.”

Slope aspect and wind persistence both show relatively low correlations with A. However, on long timescales these covariates are expected to be largely independent of temperature and precipitation, respectively, unlike density and grain size, so were included in the MLR presented here. For aspect, backscatter is minimised for a slope aspect azimuth of 11°. This has important effects along the EAIS divide: the observed A discontinuity across the divide is not present in any other covariate. Wind persistence was chosen as a covariate because of the higher correlation with A than wind speed (r =  + 0.21 vs r = −0.16, respectively), and the possibility of wind modifying the bulk backscatter via either (isotropic) surface roughening or surface grain size modification (Linow and others, 2012). Furthermore, wind persistence was chosen in order to maintain a common set of covariates for both the A and M ₂ MLRs (described below; persistence is important in that case).

Upper 1 m mean grain radius shows a reasonable positive correlation with A, as predicted by radiative transfer models (Partington and others, 2003). Upper 1 m mean density also shows a reasonable correlation with A (r = −0.54). Spatially, the A recreation based only on density closely resembles that of precipitation, as expected (higher precipitation gives lower density). Given this close relationship, there is likely little independent information added by using both precipitation and density. Indeed, performing the MLR without the grain size and density covariates results in a slight increase to the RMS residual (from 1.64 to 1.73 dB, or from 8.2 to 8.7% of the range of the A parameter).

The spatial distribution of the A parameter residual reveals a systematic A underestimation in areas of exposed rock (shown as red regions in Fig. 3c), e.g. along the Transantarctic Mountains dividing West and East Antarctica, Fimbulheimen in Dronning Maud Land, the Prince Charles Mountains lying south-west of the Amery Ice Shelf, the Shackleton Range upstream of the Filchner Ice Shelf and particularly the Ellsworth Mountains (~78.5°S, 85°W). These locations are indicated in Figure 2. This is not unexpected, since the backscatter of rock is generally higher than the surrounding polar ice (Fig. 4a in Fraser and others (2014), particularly around the Crary Mountains (~76.75°S, 117.6°W) and Mt Takahe (~112°W; 76.25°S)), and is otherwise unaccounted for in the MLR process. Other than regions of exposed rock, A appears to be underestimated in glaze-prevalent regions. This will be discussed next with the MLR process for the M ₂ parameter.

Fig. 4. MLR of covariates slope aspect (e, h), precipitation (f, i), temperature (g, j), wind persistence (k, n), upper 1 m mean grain size (l and o) and upper 1 m mean density (m, p) to recreate M ₂. (a) The long-term mean M ₂ parameter. (b) The MLR recreation of M ₂ using the aforementioned covariates. (c) The residual (RMS residual = 0.86 dB, or ~17.4% of the M ₂ range).

The MLR recreation of M ₂ is of comparatively limited success (Figs 4a and b). The recreated M ₂ parameter residual error (Fig. 4c) is 0.79 dB, or ~15.8% of the ~5 dB range. The main factor determining this poor recreation is the relatively low correlation between all covariates and M ₂ (range: |r| = 0.04 to 0.42; Figs 4k–p). Upper 1 m grain radius shows the highest correlation with M ₂ (r = −0.42). This is probably because of the enhanced near-surface grain growth in glaze-prevalent regions, i.e. regions characterised by smooth surfaces (Scambos and others, 2012). That is, grain size does not drive M ₂, but both are strongly dependent on the surface condition. Aspect (transformed for linearity, as previously described for the A parameter MLR) shows a similar high correlation coefficient, r = −0.41, followed by precipitation (transformed, as with the A MLR; r = −0.33), upper 1 m mean density (r = 0.27), wind persistence (r = 0.24) then surface temperature (r = 0.04). Slopes with a northerly aspect show a wide range of M ₂, but those with more southerly aspect are limited to low M ₂ values. Unexpectedly, the relationship between M ₂ and slope aspect is stronger than that of M ₂ and wind persistence. Regions with high precipitation (i.e. low transformed precipitation in Fig. 4j) display M ₂ values in a wide range of 0–5 dB, whereas low precipitation regions show only low M ₂ values, implying that the latter are glaze-prevalent regions as detailed by Scambos and others (2012).

Analysis of the spatial distribution of M ₂ residual (Fig. 4c) reveals that positive M ₂ anomalies occur mainly in glaze-prevalent regions (Scambos and others, 2012). Das and others (2013) determined that surface scour occurs in regions with both a low ratio of accumulation to wind speed (A/W), and a high mean slope in the wind direction (MSWD). The Das and others (2013) scour criterion is recreated here (using adjusted thresholds of A/W≤11.0 and MSWD≥0.0013 to account for the lower spatial resolution used for the calculation here), and presented in Figure 4d. The Das and others (2013) criterion is a highly nonlinear process, i.e. surface scour only appears when both A/W and MSWD exceed certain thresholds. Due to this complex nonlinearity, properly accounting for surface scour is impossible in a simple MLR analysis. Nevertheless, this result indicates that accumulation and wind processes that generate glaze-prevalent regions are important driving factors for M ₂.

In summary, the long-term mean ASCAT A parameter is successfully recreated using MLR, mainly because of the strong correlation between A and long-term mean precipitation. The same technique for recreation of M ₂ suffers due to complex nonlinearities involving surface scour/glaze in regions of high MSWD and low A/W ratio (Das and others, 2013).

It is reasonable to expect that all correlations will increase when using longer time series of atmospheric driver datasets. However, when the AMPS precipitation was substituted with the RACMO2.1 long-term mean (1979–2010) SMB dataset (Lenaerts and others, 2012), the A parameter MLR residual remained essentially unchanged, at 1.65 dB. This provides independent validation of the high quality of the AMPS precipitation dataset.

3.2. Interannual trends

It is desirable to investigate the interannual trends in the scatterometer parameters to provide information on their drivers. Alternatively, apparent interannual trends can be manifestations of the systematic burial of snow with anomalous properties (e.g. burial of a snow melt layer). To investigate the latter phenomenon, we produce a first-order estimation of time to burial to a depth of 20 m (i.e. the C band effective penetration depth). Using a Marie Byrd Land density profile (Ligtenberg and others, 2011) and the positive SMB provided by the AMPS precipitation field, a dry snow zone mean burial time of ~150 years is obtained. This is much longer than the ASCAT time series (4.5 years data used here). Thus, it is difficult to distinguish between genuine trends and burial of transient anomalous layers. In reality, subsurface snow metamorphism will reduce the contrast in properties between the anomalous layer and the typical snow properties (Linow and others, 2012), so the burial time is an upper limit for the timescale of the persistence of anomalous signals. According to observations shown in Figure 5, the timescale for backscatter recovery to near-quiescent values appears to be on the order of 1 a. Figure 5 shows an example of a perturbation of A and M ₂ from Marie Byrd Land, West Antarctica (117.5°W; 83.5°S) in 2011, likely generated by the near-concurrent (austral autumn) precipitation event shown in the lower panel (no other surface parameters showed anomalous events at this time). Over the remainder of the year, both scatterometer parameters trend toward the 2008–10 mean value, which may represent the long-term mean. There is also evidence for a similar precipitation anomaly prior to the start of the time series. Given that the apparent timescale for parameter recovery is on the order of 20% of the time series length, the possibility for ‘contamination’ of apparent trends by such events is strong.

Fig. 5. Time series of A, M ₂ and precipitation from a pixel located at 117.5°W; 83.5°S (Marie Byrd Land, West Antarctica).

Nevertheless, a significance analysis (Figs 6a and b) shows that ASCAT A and M both have extensive regions of statistically-significant trends (p-value <0.05). The same analysis for the surface datasets (Figs 6c–h) produces mixed results. Temperature, wind speed and wind persistence all show pockets of statistically-significant trends. However precipitation, acknowledged as one of the largest contributors to A on multi-annual timescales (Bingham and Drinkwater, 2000), shows no statistically-significant trend anywhere on the continent. This is likely because of the episodic nature of precipitation, which drastically increases the standard error of the mean, i.e. reduces the t-statistic.

Fig. 6. Trends in the ASCAT (a, b), AMPS (c–f) and RACMO (g, h) parameters, across the 4.5 years time series. Regions showing statistically-significant trends are shaded. Seasonal cycles (cos(1t) and cos(2t) Fourier components) have been removed in all cases.

Although difficult to quantify, it is worth noting that the spatial pattern of the A parameter trend bears reasonable resemblance to the precipitation trend, with areas of negative precipitation trend and negative A parameter trend observed in West Antarctica, Wilkes Land and western Dronning Maud Land. These are interspersed with positive precipitation and A parameter trends in eastern Dronning Maud Land/Enderby Land, and George V Land. This points to precipitation being a driver of the observed A parameter trends.

Thus, while statistically-significant trends are observed in the ASCAT data, statistical significance is not achieved with the short precipitation time series. Furthermore, apparent trends in the ASCAT parameters are strongly influenced by transient events and their subsequent recovery. In light of this result, it is concluded that a statistically rigorous investigation of A and M ₂ interannual trend forcing requires a much longer time series, ideally of both ASCAT and forcing parameters.

3.3. Annual cycle investigation

Throughout the dry snow zone, both A and M ₂ display strong annual cycles. The aim of this section is to characterise these, and to elucidate likely drivers of observed seasonality.

The annual cycle is characterised by fitting a two-term Fourier series to the A and M ₂ time series for each pixel, after removing the long-term mean and linear trend. The seasonal cycle thus takes the form:

(4) $$A(t) = {k_1}\cos (2\pi t - {\phi _1}) + {k_2}\cos (2(2\pi t - {\phi _2}))$$

where t is the time (a), k ₁ and k ₂ are the coefficients of the first and second Fourier terms and ϕ ₁ and ϕ ₂ are the phase angles (radians). The cos(1t) Fourier parameters are displayed in Figure 7. k ₁ is generally larger than k ₂ everywhere. A varies largely in-phase across the entire dry snow zone, with a maximum value typically achieved in July–September (Fig. 7b). By contrast, the phase angle for M ₂ shows two distinct regimes. There is a large region in East Antarctica with a low phase angle (i.e. a peak M ₂ in January–March), whereas coastal regions show an M ₂ peak much later (Fig. 7d). The early peak region corresponds well to maps of megadunes by Fahnestock and others (2000), and glaze-prevalent regions observed by Scambos and others (2012).

Fig. 7. Fourier coefficients of the cos(1ϕ) parameter for the seasonal cycles of A (top row) and M ₂ (bottom row).

The Fourier analysis is useful for characterising the seasonal cycle. However, when considering only the two Fourier series coefficients controlling the cos(1t) Fourier component, the number of degrees of freedom is reduced to only two per pixel. In order to retain more independent data, the determination of which driving factors influence the A and M ₂ seasonal cycle is carried out instead using 90 days low-passed time series data, i.e. four independent data points per year, or 18 independent points for the ASCAT time series. This gives a two-tailed correlations significance threshold (p <0.05) of r = 0.47. As stated in Section 1, although proof of attribution is difficult without using a radiative transfer model, a lag correlation analysis at least indicates which factors could possibly contribute to the observed variability. When guided by prior research and knowledge of physical mechanisms affecting backscatter, lag correlation analysis (including both the magnitude and sign of lag) provides important insights into the processes. Furthermore, the spatial structure and coherence of the correlation shows which relationships are widespread, and which are localised.

The highest correlation between the low-passed A parameter and the driving parameters is between A and temperature (mean dry snow zone correlation is r = −0.60; Fig. 8a). This correlation is statistically-significant throughout the majority of the dry snow zone. An A parameter lag of 25 days is required to maximise the correlation.

Fig. 8. Analysis of the drivers of the seasonal cycle. Upper panels: correlation coefficient between A or M ₂ and each driving factor, as a function of lag. Solid curve: dry snow zone mean correlation coefficient. Dashed (dotted) curve: correlation coefficient mean only for positive (negative) correlation pixels. The ‘X’ on each plot shows the lag chosen for the corresponding lower panel. Lower panel: correlation map between A or M ₂ and each driving factor, at the lag indicated in the upper panel. Statistically-significant correlations are shaded.

As mentioned in the long-term A parameter analysis, the observed anticorrelation with surface temperature (i.e. higher temperature gives lower backscatter) implies that temperature drives the A parameter via modification of the dielectric constant, rather than via grain growth processes. The lag of maximum correlation at 25 days can be interpreted as propagation of the temperature cycle vertically down through the firn. Assuming temperature propagation speed is on the order of 20 m a⁻¹ (Cuffey and Paterson, 2010), 25 days represents a propagation of temperature to a depth of ~1.4 m. This mechanism presents a plausible explanation for the observed backscatter response. Furthermore, this mechanism has been suggested in previous studies. Bingham and Drinkwater (2000) noted that the observed cycle of Antarctic backscatter was possibly explained by temperature variations, rather than precipitation cycles, since the correlation between precipitation and emissivity variability was low. This mechanism has also been explicitly investigated in model studies (Partington and others, 2003).

A also shows a strong correlation with wind speed at similar lag (A lags by 20 days, r = +0.53; Fig. 8). Wind speed also has potential links to A, mainly via (isotropic) changes to surface roughness, removal of fresh snow and reduction of surface grain size. However, given the continent-wide coherence and resemblance to the A vs temperature correlation, this relationship is interpreted to be a manifestation of the close correlation between temperature and (katabatic) wind speed, particularly on a 90 days timescale. This close correlation is also likely responsible for the reasonable correlation between A and wind persistence (Fig. 8d).

Precipitation shows very low correlation with A at all lags on this timescale (Fig. 8c). Furthermore, the structure is neither coherent across the continent, nor similar to the glaze/non-glaze structure displayed by many M ₂ correlations (details below). It is thus concluded that precipitation is not an important factor in determining the seasonal cycle of A, reinforcing the suggestion of Bingham and Drinkwater (2000).

Since both grain size and density are snow/firn properties, i.e. should be directly related to backscatter, the zero-lag correlation should reveal the most physically-realistic links between these parameters. Neither grain size nor density shows widespread correlation patterns (Fig. 8e and f), nor high mean correlation values, so the seasonal cycle is likely not driven by either upper 1 m grain size or snowpack density.

As with the A low-frequency variability analysis, M ₂ shows no strong relation to either precipitation, upper 1 m grain size or upper 1 m density (Fig. 8i, k and l). The highest correlation coefficients with M ₂ are for wind speed (r = 0.17; 30 days lag), surface temperature (r =  + 0.16; 30 days lag) and wind persistence (r =  + 0.14; 10 days lag) (Fig. 8g, h and j). As with the interpretation of the A correlations, it is important to note the strong anticorrelation between temperature and katabatic wind strength/persistence. The spatial pattern of correlation between M ₂ and these three drivers shows an area of opposite correlation in the centre of East Antarctica, as shown in the seasonal cycle characterisation in Figure 7d. This is strongly suggestive of wind speed-driven smoothing of the surface in glaze-prevalent regions, and wind-speed driven microrelief generation in accumulating regions.

As with the long-term regression analysis, the A parameter seasonal cycle is better represented by atmospheric datasets than the M ₂ parameter cycle on this timescale. Again, this is attributed to nonlinearities and complex interactions associated with wind scouring in glaze-prevalent regions.

In conclusion, the 90 days low passed A parameter dataset shows a strong correlation with both surface temperature (leading by 25 days) and wind speed (leading by 20 days). Since the temperature correlation is higher and statistically-significant over a larger area, it is likely that temperature drives A on this timescale. The high wind speed correlation is likely a consequence of katabatic wind enhancement during colder periods. The negative correlation between A and surface temperature implies that modification of the dielectric coefficient is the mechanism responsible. The M ₂ cycle is likely driven by wind speed. Glaze-prevalent regions show an anticorrelation with wind speed (wind scouring of the surface), whereas accumulating regions show a positive correlation (generation of sastrugi and microrelief).

3.4. High-frequency backscatter driver investigation

Because of the 5 day orbital sub-cycle for MetOp-A, ASCAT can generate near-complete maps of the polar regions on a 5 days basis. This is exploited here to analyse the high-frequency variability in A and M ₂, and to investigate potential drivers.

Both A and M ₂ parameters show variability at this timescale. Noise in the scatterometer parameters is expected to be much greater in this analysis than in the seasonal cycle analysis because of less temporal averaging and lower viewing angle diversity. The ASCAT radiometric RMS error of ~0.35 dB (Fraser and others, 2014) is on the same order as, or larger than, the magnitude of some high-frequency ‘events’ in the A and M ₂ parameters (e.g. Fig. 5). Insufficient azimuthal sampling, particularly south of ~78°S where only one of the two ASCAT swaths is available, may also lead to aliasing in the backscatter parameterisation process. This would most strongly manifest as noise in the high-frequency analysis. For these reasons, the correlations between A and M ₂ and the driving parameters are also expected to be lower. However, since each data point is essentially independent of the preceding and following points, DOF ≃ the number of points (n = 335 for the ASCAT time series used here). The correlation significance threshold is accordingly reduced, from r = 0.47 (in the seasonal cycle analysis where n = 18) to r = 0.11 (two-tailed, p < 0.05).

As with the seasonal cycle analysis, a lag correlation analysis is used to determine potential drivers. In order to be a potential driver, the lag between impulse and response must be either zero or positive (i.e. A or M ₂ must lag the driver). In the case of this high-frequency analysis, the lag for maximum correlation is expected to be much shorter than the seasonal cycle analysis. This analysis is shown in Figure 9. Since the density and grain size time series were available only at monthly resolution, these were excluded from the high-frequency analysis. However, these parameters are not expected to be strong drivers of the high-frequency variability, given the low correlations shown in the seasonal cycle analysis.

Fig. 9. High-frequency response analysis. Upper panels: correlation coefficient between A or M ₂ and each driving factor, as a function of lag. Solid curve: dry snow zone mean correlation coefficient. Dashed (dotted) curve: correlation coefficient mean only for positive (negative) correlation pixels. The ‘X’ on each plot shows the lag chosen for the corresponding lower panel. Lower panel: correlation map between A or M ₂ and each driving factor, at the lag indicated in the upper panel. Statistically-significant correlations are shaded.

As with the seasonal cycle analysis, the largest mean correlation between A and the atmospheric drivers is with the surface temperature parameter (mean r = −0.18; Fig. 9a). Again, the correlation is negative, implying that modification of the dielectric constant is the mechanism responsible for this relationship. At an A parameter lag of only 5 days, the penetration of the temperature anomaly signal into the snow/firn is very shallow (within ~0.27 m assuming a propagation speed of 20 m a⁻¹). This is a small fraction of the penetration depth (~20 m), however returns from the near-surface have a greater weight than deeper returns. This insufficient penetration time is likely a factor in the relatively low correlation coefficient. Nevertheless, the spatial pattern of negative correlation is widespread throughout the dry snow zone, as with the seasonal cycle analysis. Weaker positive correlations are found between the A parameter and wind speed/persistence. Again, these are considered as secondary responses to the associated temperature anomaly, and not as direct drivers of the A parameter. As with the seasonal cycle, precipitation appears to have little relationship to A at this timescale.

For the M ₂ high-frequency analysis, the wind persistence metric shows the highest correlation coefficient (r = 0.16), at zero lag (Fig. 9h). The spatial pattern shows widespread statistical-significance is achieved everywhere outside the glaze-prevalent regions on the plateau. This fits well with the conceptual model of M ₂ enhancement generated by persistent katabatic winds over a ~40 days timescale. The glaze/non-glaze discontinuity seen in the seasonal cycle analysis is not obvious here, indicating that the characteristic timescale for forming glaze-prevalent surface conditions is longer than the 90 days filter cut-off. Mean dry snow zone correlation coefficients with the other three atmospheric parameters are much lower (ranging from r = −0.05 to 0.09), indicating that wind persistence is likely the dominant driver of M ₂ variability on short timescales in non-glaze-prevalent regions.

4. DISCUSSION, IMPLICATIONS AND CONCLUSIONS

This study has documented the variability in bulk backscatter (the A parameter) and bi-sinusoidal azimuthal modulation strength (the M ₂ parameter) on the dry snow zone of the AIS. The considerable variability documented here has implications for ‘snapshot’ studies of the AIS (e.g. parameter retrievals based on a single SAR image), or studies based on sub-annual means (e.g. work based on the Seasat-A Satellite Scatterometer (Long and Drinkwater, 1994) dataset, which is only 4 months in duration; or the NASA Scatterometer (Long and Drinkwater, 1999; 9 months duration). Variability may be even larger for higher frequency instruments, since the penetration depth in dry snow is lower, meaning returns come from closer to the surface. The work presented here gives an indication of the regions in which sub-annual means and snapshots may be representative of longer-term conditions (i.e. inland areas).

The M ₂ parameter MLR analysis highlighted the need to account for the nonlinear precipitation and wind processes, which combine to form surface glaze on the East Antarctic plateau. The MLR analysis of long-term mean surface parameters and the A parameter showed that A is highly correlated with precipitation on long timescales. It is hypothesised that precipitation may be the primary driver of A, though there is some suggestion that direct modification of the dielectric constant of snow/firn via the downward propagation of temperature signals may also contribute. Owing to the high correlation between these parameters on long (>1 a) timescales, it appears that A can be used as a proxy for precipitation (or long-term mean net SMB) in the dry snow zone of the AIS (and, by extension, the dry snow zone of the Greenland ice sheet). Such studies have been conducted in Greenland (Drinkwater and others, 2001), and on small scales in Antarctica (Dierking and others, 2012), but have yet to be applied on a large scale in the Antarctic. Scatterometers such as ASCAT are well suited to this task, given their wide incidence and azimuth angle diversity, which can be used to characterise and remove anisotropy. Retrieval of long-term mean net SMB from ASCAT is the topic of future work.

Although there is potential to apply these techniques to studies of the Greenland ice sheet, recent heavy melt conditions may limit applicability. Prior to 2012, the Greenland ice sheet had a sizeable dry snow zone, however the entire ice sheet experienced surface melt in July 2012 (Nghiem and others, 2012). Thus, studies of this kind using future instruments may be impossible due to the inclusion of a strongly-scattering refrozen layer, and ASCAT prior to 2012 may be the best scatterometer with which to study the pre-2012 GrIS. The post-2012 GrIS provides an interesting way of studying the effect on backscatter of burial of a strongly-scattering layer, and ASCAT is also well suited to this task.

Statistically-significant trends were observed for both the A and M ₂ parameters. However, the observed trends potentially reflect the recovery of post-anomaly backscatter values to the long-term mean value (timescale for recovery is observed to be ~1 year). Thus, a longer ASCAT time series is desirable. While precipitation shows trends with spatial similarities to the A parameter trends, the episodic nature of precipitation time series increases the standard error, meaning significance is not achieved within a 4.5 years time series. However, given that precipitation is likely the dominant driver of the A parameter on long timescales, it stands to reason that trends in precipitation also drive trends in A. For robust trend analysis, longer time series of both scatterometer and surface datasets are required.

The seasonal cycles of A and M ₂ were characterised by fitting a two-term Fourier series to low-passed time series data. In contrast to the widespread in-phase seasonal cycle in the A parameter with a backscatter peak in July–August, the M ₂ cycle in accumulating regions is out-of-phase with the M ₂ cycle in glaze-prevalent regions where accumulating faces are interspersed with extensive regions of near-zero net SMB. In glaze-prevalent regions, the M ₂ parameter maximum is observed in January or February, i.e. out-of-phase compared with accumulating regions. This shows the large difference between surface characteristics in the two regimes: Glaze prevalent regions are smoothed or scoured by stronger winter katabatic winds, while fresh snow in accumulating regions is (highly anisotropically) roughened.

Analysis of dry snow zone mean scatterometer and surface parameters revealed the complex relationship between the seasonal cycle of A and surface temperature. In agreement with theory (Partington and Flach, 2003) and past studies (Bingham and Drinkwater, 2000), it appears that the annual cycle of A is primarily driven by the vertical propagation of the surface snow/firn temperature cycle, with a maximum A typically occurring in July or August, lagging the temperature cycle by ~25 days.

It is generally held that subsurface grain size (within the penetration depth) is an important driver of microwave backscatter (Ulaby and others, 1996; Partington and Flach, 2003). While this appears to hold to some extent on a multi-year timescale (upper 1 m grain size shows a reasonable correlation with the long-term mean A parameter (r =  + 0.46)), it appears that precipitation is much more important (r =  + 0.78, after transformation). Furthermore, this study has shown that the seasonal cycle of grain size variability shows no widespread correlation with the A parameter cycle. Instead, modification of the dielectric properties of the snow/firn by temperature cycle propagation dominate the A parameter seasonal cycle.

Analysis of 90 days high-passed time series again reaffirmed the strong links between the A parameter and temperature anomalies. For the M ₂ parameter, as with the seasonal cycle analysis, wind parameters (especially persistence) showed the highest correlation, particularly in non-glaze-prevalent regions, as expected.

Polar firn displays much lower backscatter variability than the return from the ocean, especially on short timescales. Thus, this dataset provides an ideal means for validation of the newly-commissioned MetOp-B ASCAT, and future ASCAT instruments. Furthermore, the launch of identical ASCAT instruments is important for mission continuation, and will likely provide the first decade-long scatterometer dataset. Studies of the ice sheets will benefit greatly from this longer time series, so ASCAT is set to become an important instrument for ice sheet change detection and attribution.

Given the difficulties associated with change attribution, this type of study would benefit from the inclusion of scatterometer data of different frequency and polarisation, provided the ability to characterise and remove anisotropy is on par with ASCAT (i.e. azimuth and incidence diversity is maintained). Higher frequency instruments have much lower penetration depths, which should enable better attribution of surface changes. Inclusion of nadir active microwave instrument information (e.g. radar altimeters) can also contribute independent information which will be helpful for attribution. The conclusions drawn from this study may also impact on the interpretation of radar altimetry. Thus, a broadening of this kind of study to include other instruments is a priority.