2.2.1 Calculating RMS height, R

As noted, subglacial-roughness information can be determined via the statistical analysis of vertical variation in along-track bed topography (e.g. Siegert et al., 2004, 2005; Taylor et al., 2004; Rippin et al., 2006, 2011, 2014; Bingham and Siegert, 2007, 2009; Bingham et al., 2007, 2017; Li et al., 2010; Rippin, 2013). The most prevalent method in glaciological literature employs spectral methods to do this (i.e. the application of fast Fourier transforms – FFTs – first employed in glaciology by Hubbard et al., 2000, and for the Antarctic Ice Sheet by Taylor et al., 2004). Alternative space-domain methods exist, however, and are frequently used within earth and planetary sciences (Shepard et al., 2001; Smith, 2014).

Here, the first metric for subglacial roughness we present, topographic roughness (or R), is quantified by the RMS height in along-track topography (RES sampled bed elevation). The RMS height (referred to also as the standard deviation of bed elevation; e.g. Rippin et al., 2006, 2014) provides several benefits over the use of FFTs. First, it enables the collation of all CReSIS survey campaigns despite variable sample spacing (horizontal resolution) without requiring along-track interpolation or re-sampling of data. Second, this method allows the use of a shorter length scale than FFTs, not only facilitating subsequent anisotropic analysis at crossovers (Sect. 2.2.2) but also providing a finer-scale roughness information. A final advantage is that RMS height calculations are unit-preserving (i.e. quantifying variation at the bed in units of metres), providing a more physically intuitive metric. More critically, however, the spatial distribution of roughness values quantified by FFT and RMS height methods have been noted to be similar (Rippin et al., 2014; Falcini et al., 2018).

Sampled bed and surface elevations were obtained from all the available CReSIS RES surveys between 1993 and 2016. Where applicable, data were filtered using the provided quality flags denoting the confidence of the bed pick accuracy (Paden, 2017), ensuring that only bed elevations with “high” confidence were used; however, as RES data obtained during OIB campaigns prior to 2008, with the exception of the reprocessed “2006 TO” survey, do not include quality flags, all available sampled bed elevations were used. More-recent surveys (post-2006) have an increased along-track sampling resolution (approximately twice that of previous campaigns), owing to SAR processing and multi-looking stages in data preparation. This results in an “overlap” between consecutive samples of bed elevation (John Paden, personal communication, 2016). Therefore, to ensure that only independent measures of bed elevation are used for roughness calculation, data from these campaigns are rarefied (to include every other sample point).

Topographic roughness, R, is given by

where n is the number of sample points, z(x_i) is the height of the surface point at point x_i, and $\stackrel{\mathrm{‾}}{z}$ is the mean height of the profile over all x_i values. RMS height R was calculated using a window length or bin size, L, of 200 m, using all recorded bed elevations, regardless of spatial density within the bin, provided that n≥3. R is given for the spatial midpoint of each window. Regions of greater roughness, quantified by a larger variation in bed elevation within the window, have greater R values. An example of R calculated along track using the sampled bed elevation is presented in Fig. 2. It should be noted that not all bins have a constant n due to the variation in sampling regime, the resolution of the different radar instruments, and data quality. R was not calculated for bins where n<3. Although n is small (mean is ∼8), we obtain a large sample size for calculating roughness statistics through the repeated sampling over multiple bins (using repeat flight tracks).

Figure 2Along-track example of calculated topographic roughness (R). This demonstrates the length scale (200 m) over which R is calculated from sampled bed elevation. Grey bars depict high values of R associated with subglacial step changes in elevation (“cliffs”); the limitation of interpreting topographic roughness in these regimes is discussed in Sect. 4.4.

Download

L=200 m was chosen to enable the finest scale of R to be quantified whilst maintaining the largest spatial coverage of the resultant metric by using all available survey data. It is possible to quantify R at a finer scale using only more-recent survey data; however this is at the expense of reduced spatial coverage for a length scale not less than 100 m (limited by the along-track sample spacing). Changes in L influence the quantification of roughness as a result of the self-affine (fractal) scaling behaviour of subglacial terrain: as L increases, bed profiles with a steeper slope tend to become more rough relative to those with shallower slope (see Figs. 3, 1, and 2a, respectively, in Shepard et al., 1995, Shepard and Campbell, 1999, and Jordan et al., 2017).

2.2.2 Filtering R with respect to ice surface velocity

To evaluate, and more completely understand, how the spatial distribution of subglacial roughness influences or perhaps is influenced by ice-sheet motion, we compare R to ice surface velocities. We use the InSAR-derived MEaSUREs velocity mosaic (Joughin et al., 2016, 2017) over the entire GrIS. This mosaic helps with capturing long-term information (using 1995–2015 observations) regarding flow configuration, minimising inter- and intra-annual variation in both ice speed and direction. As R is quantified using 2 decades of RES data, we assume an inherent constancy in roughness over time. The MEaSUREs data provide magnitude ($\left|v\right|$; speed) and direction at a 250 m resolution (Fig. 3a and b); however, for our analysis, we performed a bilinear aggregation (to 1000 m) in order to smooth small-scale variation or noise.

Figure 3Observed surface ice velocity characteristics of the GrIS used in the filtering of R. (a) InSAR-derived ice surface speed (velocity magnitude; m a⁻¹; Joughin et al., 2016); regions of fast ($\left|v\right|>\mathrm{50}$ ma⁻¹) and slow ($\left|v\right|<\mathrm{5}$ ma⁻¹) flow are demarcated by the black and white contour lines, respectively. (b) Flow direction of ice surface from (a); coloured pinwheel denotes direction of surface ice flow, where north is at the top of the page. (c) Radar sounding surveys as in Fig. 1, filtered for alignment with surface flow direction (b); flight tracks are categorised as aligned either parallel (R_∥) or perpendicular (R_⟂) to surface flow direction (with a ±20^∘ threshold) for the analysis of topographic-roughness anisotropy.

With regard to ice surface flow speed, we delineate regions of “fast” ($\left|v\right|\ge$ 50 ma⁻¹) and “slow” ($\left|v\right|\le \mathrm{5}$ ma⁻¹) flow (Fig. 3a). In regions where $\left|v\right|$ exceeds 50 ma⁻¹, ice is likely to be decoupled from the bed (i.e. sliding), as this speed cannot be achieved by internal deformation alone (MacGregor et al., 2016; Stearns and van der Veen, 2018). As we have noted above, basal traction, a principal constraint on basal sliding (Weertman, 1957), may be influenced by subglacial roughness (Siegert et al., 2004, 2005; Bingham et al., 2017). Second, where ice motion is limited in slow-flowing regions, rates of basal erosion are minimal and, thus, the influence on subglacial topography is reduced (Bingham and Siegert, 2009).

As R is quantified along track, there is an inherent directionality in its characterisation of the subglacial environment. To assess anisotropy at the bed, with particular reference to ice motion, we classify R through its alignment with local flow direction. Sample windows were filtered for their linearity to remove measures of R over corners and bends (with a deviation ≥10 %; after Bingham et al., 2015) in RES flight lines. Roughness bins were then filtered by their alignment to local surface ice flow direction (Fig. 3b) with a 20^∘ threshold; Fig. 3c shows classified measures of R aligned perpendicularly R_⟂ or in parallel R_∥. From this, we draw conclusions based on the relationship between subglacial roughness and the speed of overlying ice. It should be noted that we only use contemporary ice velocity observations in this study; although flow configuration is likely to have remained largely constant, the surface speed will have changed through time. In terms of fractional error and uncertainty regarding ice speed and flow direction, and thus the classified alignment of roughness bins, slow-flowing regions represent the worst-case scenario. General error propagation formulae for independent velocity vectors (v_x and v_y) were applied, giving a mean error of 0.51 m a⁻¹ for $\left|v\right|$ and 14.55^∘ for direction, of which the latter is less than the larger angular threshold used to classify the alignment of roughness bins with respect to velocity.

Where coincident measures of R_⟂ and R_∥ are available (the near-orthogonal crossovers between flight lines), the degree of anisotropy can be calculated. This is achieved through a normalised difference ratio, herein termed “anisotropy ratio” (Smith, 2014), given by

Here, using Ω, we map the distribution of roughness anisotropy across the GrIS and assess the relationship between $\left|v\right|$ and Ω in both fast- and slow-flowing regions. Values of Ω are interpreted such that −1 dictates a complete dominance of smoothness parallel to flow direction (perhaps as a result of flow-aligned features), where +1 is a dominance of smoothness perpendicular to flow (i.e. parallel roughness) and values of ∼0 indicate roughness isotropy.

2.3.1 The abruptness (peakiness) of the bed-echo waveform

Bed-echo waveform properties are related to electromagnetic scattering from the glacier bed and, hence, also provide information about subglacial roughness (Oswald and Gogineni, 2008; Oswald et al., 2018; Jordan et al., 2017). Radar bed-echoes range from sharp pulse-like returns (associated with specular reflections from a smooth glacier bed) to echoes that have a trailing edge that extends greatly over the original pulse length (associated with diffuse scattering from a rough glacier bed). A convenient way to parameterise the relative spread of the bed-echo waveform is to use the waveform abruptness parameter, defined by

where P_peak is the peak power of the bed echo and P_agg is the aggregated (integrated) power over the echo envelope (Oswald and Gogineni, 2008; Jordan et al., 2017). Three examples of bed-echo waveforms, and their abruptness values, are shown in Fig. 4c. Higher values of A are associated with specular reflections, and lower values are associated with diffuse scattering. However, the maximum value for A (A_max; which is constrained by the ratio of the image sample rate to depth-range – vertical – resolution; Jordan et al., 2017) can differ between different CReSIS field seasons, with values ranging between 0.5 and 0.8 (as determined empirically from the abruptness distribution). Since the RES bed echo results from a superposition of along-track and cross-track energy, the abruptness is a near-isotropic or isotropic parameter (Young et al., 2016) and therefore obscures information regarding the anisotropy of the glacier bed.

Figure 4Estimation of scattering-derived roughness and data combination for bed-echo peakiness. (a) Abruptness as a function of RMS height for different CReSIS field seasons: solid black curve – 2010 DC8; long dashed blue curve – 2011 TO, 2011 P3, 2012 P3, 2013 P3, and 2014 P3; dotted red curve – 2006 TO; dashed green curve – 2005 TO; and solid grey curve – 2008 TO and 2009 TO. (b) Peakiness index as a function of wavelength-scaled RMS height. (c) Example bed-echo waveforms, their abruptness (A), peakiness index (Λ), and scattering-derived roughness (ξ). The plots are for the 2011 P3 field season, which has maximum A∼0.65.

Download

The procedure used to extract the bed-echo abruptness from CReSIS Level 1B data is outlined in Jordan et al. (2017). Briefly, this consists of the following three steps. First, CReSIS Level 2 picks are used as initial estimates for the depth-range bin of bed-echo power peak. Second, a local retracker is used to locate peak power. Third, the power is integrated over the bed-echo envelope, applying a “quality control” measure such that the peak power is 10 dB over the noise floor. This final step results in some regions, primarily in southern Greenland, having reduced coverage (see Fig. 1b in Jordan et al., 2018).

2.3.2 Estimating fine-scale roughness and the “peakiness index”

The scattering of the radar pulse at the glacier bed is underpinned by the physics of electromagnetic diffraction (Berry, 1973; Ulaby et al., 1982). As bed roughness increases, the radar pulse is scattered over a greater range of angles; this results in a decrease in peak returned power and an increase in the trailing edge of the echo. The mathematical formulation of this relationship depends on the physical model for electromagnetic interference (phase coherence or incoherence) and the statistical model for the subglacial interface (Berry, 1973; Peters et al., 2005; Haynes et al., 2018). The most commonly employed scattering model for the RES of glacier beds assumes phase-coherent interference, “smoothly undulating” Gaussian statistics for RMS roughness and radial isotropy (Berry, 1975; Peters et al., 2005; MacGregor et al., 2013; Grima et al., 2014; Schroeder et al., 2015). We employ this scattering model for two objectives: firstly as a way of estimating “fine-scale roughness” from the abruptness and secondly as a way of combining the abruptness for different radar systems to derive an (approximately) system-independent peakiness index (Λ).

Following a similar approach to that described by Schroeder et al. (2015) and Jordan et al. (2017), under assumptions of energy conservation, the scattering model can be used to predict the relationship between A and RMS height ξ (“fine-scale” roughness). In this context ξ is not strictly equivalent to the values obtained from topography, and a length-scale separation is performed with respect to a reference plane (Berry, 1973). The relationship between A and ξ is given by

where

denotes the RMS phase variation, with A_max being the maximum abruptness, I₀ being a 0th-order Bessel function of the first kind, f_c being the centre frequency of the radar pulse, c being the vacuum speed of the radar pulse, and ϵ_ice=3.15 being the relative dielectric permittivity of glacier ice (Peters et al., 2005). Since the radar wavelength in ice is ${\mathit{\lambda }}_{\mathrm{ice}}=c/f_{\mathrm{c}}{\sqrt{\mathit{ϵ}}}_{\mathrm{ice}}$, Eq. (5) can be expressed as

and hence ξ is scaled by the radar wavelength in ice (either 0.87 or 1.13 m for the 195 and 150 MHz systems, respectively). There are therefore 2 degrees of freedom in Eq. (4) that can vary for different CReSIS field seasons: A_max and f_c. The different parameter combinations are shown in Fig. 4a, and from these relationships it is possible to estimate ξ from A (and thus obtain a measure of fine-scale roughness that is similar between different radar systems). However, since the values of f_c and A_max differ between field seasons, a crossover bias is present for “raw” abruptness values. In order to combine abruptness data, we back-substituted the value of ξ to obtain the value of A as if it were the most spatially extensive radar system (the blue curve in Fig. 4a) and then linearly rescaled amplitude on the interval [0,1] (using A∕A_max) to give the peakiness index (hereafter referred to as Λ). These steps combine the measurements via the system-independent relationship that is modelled between Λ and wavelength-scaled RMS height ξ∕λ (Fig. 4b).

The inter-season data combination was validated by performing crossover analysis for ξ and Λ, with the allowed tolerance for the crossover bias set to 5 % of the parameter range. RES data that do not meet this criterion (primarily the older ICORDS data but also the data from 2010 P3 season, which is known to have noise-floor issues; Paden, 2017) were discounted completely from analysis. Although the data combination scheme employed here, across CReSIS platforms, is seen to work well, it should be noted that combining data from multiple instruments, particularly those with a large difference in centre frequencies, may not be so effective.

It is important to note that obtaining ξ from Eq. (4) is just one way of estimating fine-scale roughness. Self-affine (fractal) statistics (Shepard and Campbell, 1999) can also be applied to scattering models of glacier beds (as in Jordan et al., 2017). Additionally, in reality, fine-scale roughness is anisotropic, as revealed by the “specularity” scattering metric (Schroeder et al., 2013, 2014, 2015; Young et al., 2016). We therefore recommend that ξ should be interpreted in a qualitative manner, with lower values indicating “fine-scale smooth” and higher values indicating “fine-scale rough” regions of the glacier bed. In regions of complex bed topography, and in particular at outlet glacier regions, off-nadir scattering may adversely influence the signal and lead to a breakdown in the interpretation of this metric (see Sect. 4.4). Fine-scale roughness that relates to radar scattering can also be estimated from the statistical distribution in peak bed-echo power (Neal, 1982; Grima et al., 2014).

3.1.1 Topographic roughness, R

Across the ice sheet, unfiltered (with respect to surface flow direction) R shows clear spatial heterogeneity (Fig. 5a); coherent signals, representing contiguous regions of both “smooth” (low R values) and “rough” (high R values) beds, are visible. Generally, the margins of the ice sheet contain the roughest beds, whereas the interior is notably smooth. Ice-sheet-wide, the lowest values of R are observed in the north and north-west of the island. However, localised to the main “trunks” of Petermann and Humboldt glaciers (PG and HG, respectively, in Fig. 5a), at the point of highest $\left|v\right|$ immediately before the grounding line, small patches of smooth bed are observed. Broadly speaking, across the ice sheet, fast-flowing regions exhibit rough beds; though, as exemplified in the north and north-west, this behaviour is somewhat spatially variable at the perimeter of Greenland. Notable examples of contiguous smooth beds near the margins include the following: north-west of the Camp Century (CC) drilling site, in the vicinity of Ìngia Isbræ (II; north of Rink Isbræ), and a region near the outlet of the North East Greenland Ice Stream (NEGIS; as marked in Fig. 5a). The highest values of R trace the Caledonian fold belt mountain range (formed in ∼420 Ma ago; Henriksen, 2008) and the deep inland fjord-like systems along the eastern and south-eastern margins of the island (Fig. 5a).

Figure 5Topographic roughness (R) across the GrIS. (a) R unfiltered by flow direction. (b) R_⟂. (c) R_∥. (d) and (e) show spatial interpolation of (b) and (c) to a width of 20 km, respectively, for improved visualisation. Locations for Ìngia Isbræ (II), NEGIS, Petermann Glacier (PG), Humboldt Glacier (HG), and Camp Century (CC) drilling site, as referred to in the text, are marked.

Figure 5b and c present directionally filtered values for topographic roughness, aligned perpendicular (R_⟂) and parallel (R_∥) to ice surface flow direction, respectively. For improved visualisation (and visual analysis) only, maps for R_⟂ and R_∥ were interpolated (using inverse-distance weighting) to a limit of 10 km (Fig. 5d and e). This interpolation distance is representative of the average track spacing used in the “gridded” airborne sampling regimes in fast-flowing regions (e.g. surrounding Jakobshavn Isbræ – JI – and Petermann Glacier – PG; see Fig. 1). Initial comparison shows a marked difference between R_⟂ and R_∥, most notably in fast-flowing regions. Across the ice sheet, the bed is observed to be smoother parallel to flow. In the ice-sheet interior (where $\left|v\right|<\mathrm{50}$ ma⁻¹) the subglacial environment is mostly smooth in both directions (i.e. isotropic). However, in the south of the ice sheet we observe more distinct differences between R_⟂ and R_∥ values (see Sect. 3.2). Overall, R_∥ exhibits more uniform roughness values across fast- and slow-flowing regions, particularly within the north and west, whereas R_⟂ presents a notable difference between fast- (rough) and slow-flowing regions (smooth).

We observe a similar spatial distribution of unfiltered R (Fig. 5a) to those previously quantified for Greenland using an RMS residual technique (see Fig. 4 in Layberry and Bamber, 2001) and through a frequency-domain approach (FFTs) undertaken at a much larger length scale (3200 km; see Fig. 1 in Rippin, 2013); in these studies, general conclusions for a smooth interior and rough margin were made. Rippin (2013) additionally notes a localised smooth bed underlying the trunk of Petermann Glacier, whereas Layberry and Bamber (2001) note a smooth basin for both Humboldt and Petermann glaciers. However, as these studies do not filter with respect to surface flow direction, they do not reveal roughness anisotropy in R across the ice sheet.

3.1.2 Scattering-derived roughness, ξ

Figure 6a presents the spatial distribution of scattering-derived subglacial roughness, ξ, for the GrIS. As noted in Sect. 2.3, these values are inversely correlated to Λ (Fig. 6b) due to the scattering model relationship. The spatial distributions observed within scattering-derived roughness are broadly similar to those observed for unfiltered R, including a notable link between fast flow and high values of ξ (rougher beds). Regions that present the smoothest subglacial environments also reflect those mentioned above, notably the vicinity of the CC drilling site, a coherent patch south-east of Petermann Glacier, towards the outlet of the NEGIS, and at Ìngia Isbræ (marked in Fig. 6a). Low (smooth) values of ξ are also observed along the central ice divide. In contrast to measures of R, however, and concordant with the broad-scale relationship of ξ to $\left|v\right|$, the fastest-flowing trunks of Humboldt and Petermann glaciers contain rougher beds (HG and PG, respectively, in Fig. 6a). Other differences between topographic and scattering-derived roughness include a corridor of high ξ extending south of Petermann Glacier and across the ice divide (see Fig. 6a) as well as a generally more “mixed” roughness behaviour in the ice-sheet interior.

Figure 6Scattering-derived roughness (ξ) across the GrIS. (a) ξ. (b) Non-dimensional Λ (peakiness index) determined from the bed-echo waveform. Locations for Ìngia Isbræ (II), NEGIS, Petermann Glacier (PG), Humboldt Glacier (HG), and Camp Century (CC) drilling site, as referred to in the text, are marked.

3.2.1 Ice-sheet scale

Owing to the isotropic nature of ξ, we limit more comprehensive assessment of the relationship between contemporary ice velocity and subglacial roughness to topographic roughness, R, only (undertaken using all calculated R bins). Figure 7 presents an assessment of the relationship between R with respect to surface ice flow direction in fast-flowing regions ($\left|v\right|>\mathrm{50}$ ma⁻¹). The difference in distributions between R_⟂ and R_∥ (Fig. 7a and b) indicates that roughness perpendicular to flow direction is greater (i.e. more rough; mean is 9.39 m compared to 6.27 m) and exhibits higher variance (92.21 m² compared to 43.02 m²).

Figure 7Relationship between R_⟂ and R_∥ and surface ice velocity for fast-flowing ($\left|v\right|>$ 50 ma⁻¹) regions of the GrIS. (a) and (b) present distributions R_⟂ and R_∥, respectively. (c) and (d) show mean ice surface speed, $\left|v\right|$, calculations for logarithmic R bins (at 0.25 m intervals). This is a linear–log plot, where the limit of the horizontal axis (R) is 10^1.25 m, noted by the dashed black lines in (a) and (b). It should be noted that the vertical exaggeration of these two plots is constant. Mean values for ice surface speed are used here to facilitate direct comparability to previous work presented in Lindbäck and Pettersson (2015). Colours here are consistent with Fig. 3c for alignment with surface flow direction.

Download

Calculated mean ice surface speeds ($\left|\stackrel{\mathrm{‾}}{v}\right|$) for logarithmic bins (at 0.25 intervals) of R_⟂ and R_∥ are shown in Fig. 7c and d, respectively. A marked difference between the calculated ice speed averages is observed. For all bins of R_⟂, $\left|\stackrel{\mathrm{‾}}{v}\right|$ is seen not to exceed 250 ma⁻¹, whereas the lower bound for $\left|\stackrel{\mathrm{‾}}{v}\right|$, calculated for R_∥, is >350 ma⁻¹. This is most likely a result of a greater spread in values of $\left|v\right|$ for parallel roughness bins; however, it is notable that this scaling relationship is broadly in agreement with those previously observed in regional studies in Antarctica (Bingham and Siegert, 2007) and Greenland (Lindbäck and Pettersson, 2015). Additionally, if we are to assume that $\left|v\right|$ increases towards the glacier terminus (or grounding line), the exhibited scaling relationship for R_∥ is in agreement with previous studies where roughness is observed to decrease (Bingham and Siegert, 2007, 2009). Increasingly smooth beds parallel to flow direction, therefore, are indicative of enhanced ice surface speed. The limit to which this relationship holds is R=10^1.25 (also delineated in distribution histograms by the dashed black line in Fig. 7a and b). This value is the approximate upper limit of R that can reasonably quantified using Eq. (1) (Sect. 4.4). Conversely, a weak positive relationship is observed between R_⟂ and mean ice surface speed (Fig. 7c). R_∥, however, exhibits a strong negative exponential scaling relationship with mean ice surface speed (Fig. 7d), which is statistically significant above the p=0.001 confidence level.

Figure 8Calculated anisotropy ratio for R. (a) Anisotropy ratio, Ω, where values of −1 dictate a dominance of smoothness parallel to flow direction, values of +1 indicate a dominance of smoothness perpendicular to flow (i.e., parallel roughness), and values of 0 indicate isotropy. (b) and (c) present mean ice surface speed, $\left|v\right|$, calculated for anisotropy ratio bins (at 0.1 intervals) for slow- and fast-flowing regions, respectively.

Figure 8a presents the spatial relationship of the anisotropy ratio (Ω) across the ice sheet, where coincident values of R_⟂ and R_∥ are quantified. It is clear that fast-flowing outlet regions (the ice-sheet margins) are generally more smooth parallel to surface flow direction (where $\mathrm{\Omega }\to -\mathrm{1}$). In the ice-sheet interior a more varied, or random, distribution in Ω is apparent. Mean ice surface speed for bins of Ω at 0.1 intervals, in fast- and slow-flowing regions (Fig. 8b and c), reinforces this observed spatial relationship in subglacial roughness. A strong linear relationship with regards to $\left|v\right|$ is exhibited within fast-flowing regions, whereas in regions of slow flow no such relationship is observed.

3.2.2 Fast-flow regions and outlet glaciers

To assess any spatial heterogeneity in the exponential scaling relationship between ice flow and R_∥, local regions of fast flow were selected for closer analysis. These regions are centred around major outlet glaciers (Fig. 9) and, where possible, encompass only individual outlet glaciers (e.g. Humboldt – Region 1, Petermann – Region 2, and Kangerlussuaq – Region 4); however, where outlet glaciers are in close proximity, wider regions of fast flow were assessed (i.e. Jakobshavn+ – Region 6). Regionally, we observe the same exponential scaling relationship as exhibited ice-sheet-wide (see Fig. 10). The calculated regression line for each region is statistically significant at, or above, the p=0.01 confidence level, with the exception of Region 3 (encompassing NEGIS), at p=0.05. A marked difference in the regression gradients is also observed, spanning 4 orders of magnitude: Region 3 exhibits the shallowest gradient ($-\mathrm{1.01}\times {\mathrm{10}}^{-\mathrm{1}}$ a⁻¹), and Regions 4 and 5 exhibit the steepest gradient ($-\mathrm{9.39}\times {\mathrm{10}}^{-\mathrm{4}}$ a⁻¹ and $-\mathrm{7.66}\times {\mathrm{10}}^{-\mathrm{4}}$ a⁻¹, respectively). Echoed by the shallow regression gradient and the lower confidence level of statistical significance, the NEGIS (Region 3) also exhibits the lowest r² value (0.35). As previously described, both unfiltered R and ξ values reveal a contiguous smooth-bed signal, aligned near-perpendicular to flow direction (marked in Figs. 5a and 6a; further described in Sect. 4.1 and 4.3). Downstream from this, a coincident increase in subglacial roughness and $\left|v\right|$ is observed. Additionally, there is a notable sampling bias in the radar sounding across Region 3, where fewer tracks are aligned parallel to the flow direction (Figs. 3c and 9b). Together, these factors may be responsible for the weaker scaling relationship observed here between $\left|\stackrel{\mathrm{‾}}{v}\right|$ and R.

Figure 9Local subsets of R_⟂ and R_∥ in fast-flowing outlet glacier regions. Interpolated R_⟂ and R_∥ (as Fig. 5) are shown for the fast-flowing regions of (a) Humboldt (1) and Petermann (2) glaciers, (b) the North East Greenland Ice Stream (NEGIS; 3), (c) the north-western (NW; 7) fast-flow region, (d) Kangerlussuaq (4), (e) Jakobshavn Isbræ and surrounding glaciers (6), and (f) Helheim (5). The location of these regions is inset, where regions are colour-coded for further analysis (see Fig. 10).

More interestingly, two distinct groups are observed, showing a clear separation in regression slope gradients (Fig. 10). The first group (see bottom; Fig. 10) is mostly homogenous in terms of its regression slopes (i.e. the relationship between roughness and ice surface speed here is broadly similar). However, Region 4 and 5 in south-eastern Greenland (Kangerlussuaq and Helheim, respectively) exhibit a marked increase in gradient, indicative of a stronger scaling relationship at these sites.

Figure 10Relationship between R_∥ and ice surface speed $\left|v\right|$ for fast-flowing outlet glacier regions. This is a linear–log plot as per Fig. 7d, depicting the calculated mean ice surface speed, $\left|v\right|$, for logarithmic R bins (at 0.25 m intervals) at each of the seven regions shown in Fig. 9; the gradient of the linear model (with units of a⁻¹) for each region is shown in ascending order (less negative).

4.3.1 Camp Century

Figure 11 presents one such contiguous region of smooth beds in the vicinity of the CC drilling site: where an increase in ξ is observed towards the east and south-east, near Humboldt Glacier. Fast-flowing regions have been masked, owing to the isotropic nature of scattering-derived roughness and the anisotropic behaviour of topographic roughness outlined above (Sect. 3.2). As the bed is likely frozen in this region (MacGregor et al., 2016), where we also observe a high elevation plateau and slow-flowing ice (and a local ice divide), it is not feasible to interpret this signal as simply the presence of electrically deep basal water. From the knowledge, albeit limited, of subglacial geology in this region (see Fig. 1 in Dawes, 2009), we propose that this signal (of low ξ) is in fact caused by a non-deformable bed, related to underlying geology on which there is little to no sediment. This bed, reminiscent of preglacial erosion surfaces observed in Antarctica (Rose et al., 2015), is also likely to have been largely untouched by long-term glacial erosion.

Also observed in this region are elevated ξ values coincident with the Hiawatha impact crater (Kjær et al., 2018), associated with channelised features (triangle; Fig. 11). Whilst higher values of ξ may well be due to the interference from off-nadir echoes (as explained above; see Sect. 4.4), it is plausible that, by contrast, this may be a marker for a soft bed (i.e. presence of deformable sediment), as a result of enhanced sediment transport.

4.3.2 Igneous intrusion, Petermann Glacier

Figure 12 depicts scattering-derived roughness and bed elevation near Petermann Glacier, north-west Greenland. East of the streaming ice and bounded to the north and east by the palaeo-fluvial “mega-canyon” (Bamber et al., 2013), we observe a contiguous low-ξ region where surface flow speed is <50 ma⁻¹. This signal is observed coincident with a local topographic high (with a prominence of 300 m in elevation), which, unlike the surrounding topography, is largely left unmarked or dissected by bed channels. Previous geophysical interpretation, using both gravity and magnetic anomalies derived from OIB data (see Fig. 2 in Tinto et al., 2015), has established this unit as an intruded igneous body. The unaltered nature, and geological interpretation, of this feature further lend credibility to our interpretation of low ξ values as denoting a hard bed. Additionally, recent assessment of the basal thermal state, and basal water prediction derived from RES, suggest that this region is not predominantly “wet” (MacGregor et al., 2016; Jordan et al., 2018; Chu et al., 2018), further indicating that the interpretation of water ponding is unlikely to hold here.

4.3.3 Volcanic province, central Greenland

Well-constrained by exposed geology at the ice-free margins of Greenland (bounded west and east by Ìngia Isbræ and Geikie Plateau, respectively) is the presence of a volcanic province from the Palaeogene (Fig. 13; see also Fig. 1 in Dawes, 2009) under the inland ice in central Greenland. However, the exact extent of the presence of the underlying basaltic rocks cannot be accurately determined (Dawes, 2009). At each margin of the GrIS where $\left|v\right|$ is <50 ma⁻¹, we see good spatial agreement between ξ and the mapped volcanic province. If we are to conclude that low values of ξ delineate a hard bed, it may be possible to redraw the boundary of the volcanic province further inland from the western margin (Fig. 13). The eastern end of this “smooth” region is spatially correlated with elevated levels of geothermal heat as a result of the long-term tracking of the Icelandic hotspot, a relatively thin lithosphere, and an underplated body, discussed by Rogozhina et al. (2016) and Martos et al. (2018).

4.3.4 Delineating deformable and non-deformable beds

In many assessments of subglacial roughness in Antarctica, smooth beds have been associated with the presence of weak sediment layers beneath fast-flowing outlet glaciers (e.g. Bingham and Siegert, 2009; Rippin et al., 2014; Bingham et al., 2017) (Sect. 4.2), whereby the deformation of this sediment does not only exert important spatial controls upon the onset of fast flow but also upon ice speed in Antarctica (Alley et al., 1986; Peters et al., 2006; Siegert et al., 2016) and potentially in Greenland (Bougamont et al., 2014; Stearns and van der Veen, 2018). It is clear, however, from the regions previously described, that fast flow is not a necessary condition for large, coherent regions of “smooth” bed. Where basal conditions are not indicative of enhanced ice flow (i.e. slow-flowing, cold-based regions), we suggest that low values of ξ are indicative of a non-deformable, hard bed. However, this will only work well away from complex terrain (i.e. regions of low relief; see Sect. 4.4). Although we reject that such contiguous signals as are evidence of ponded basal water, or indeed basal thaw, it is plausible that small-scale patches of high abruptness values (high Λ and low ξ values) could still be interpreted this way.

If we extend our conclusion that ξ may be used to demarcate underlying hard beds, focus should then be drawn to regions where deforming basal sediment and sediment transport is likely to take place. As discussed, the majority of fast-flowing outlet regions exhibit high values of ξ (Fig. 6a; unfiltered R – Fig. 5a), which, by R_∥, we interpret as exhibiting basal “streamlining” influenced by $\left|v\right|$ (akin to that observed by Bingham et al., 2017, albeit at a different scale). This, alongside recent evaluation that states that many of Greenland's outlet glaciers may be driven by the availability of basal deforming sediment (Bougamont et al., 2014; Stearns and van der Veen, 2018), suggests that high values of ξ are a proxy to demarcate deformable beds.