Radar diagnosis of the subglacial conditions in Dronning Maud Land, East Antarctica

In order to better understand the spatial distribution of subglacial environments, ground-based radar profiling data were analyzed for a total distance of ∼ 3300 km across Dronning Maud Land, East Antarctica. The relationship between geometrically corrected bed returned power [ P c bed]dB in decibels and ice thickness H was examined. When H is smaller than a critical value that varies according to location, [P c bed]dB tends to decrease relatively smoothly with increasingH , which is explicable primarily by the cumulative effect of dielectric attenuation within the ice. However, at locations whereH is larger than the critical H values, anomalous increases and fluctuations in [ P c bed]dB were observed, regardless of the choice of radar frequency or radar-pulse width. In addition, the amplitude of the fluctuations often range 10∼ 20 dB. We argue that the anomalous increases are caused by higher bed reflectivity associated with the existence of subglacial water. We used these features to delineate frozen and temperate beds. Approximately two-thirds of the investigated area was found to have a temperate bed. The beds of the inland part of the ice sheet tend to be temperate, with the exception of subglacial high mountains. In contrast, the beds of coastal areas tend to be frozen, with the exception of fast-flowing ice on the subglacial lowland or troughs. We argue that this new analytical method can be applied to other regions.


Introduction
Subglacial environments of polar ice sheets are characterized by mass and energy transfers between the ice and its substrate of bedrock, sediment or water. Determining the distribution of water at the ice-sheet bed is crucial in many disciplines of polar science, such as the study of ice motion across tributaries of fast-flowing glaciers towards the inland (e.g. Bamber et al., 2000Bamber et al., , 2006Bell, 2008;Bennett, 2003;Pattyn et al., 2005;Rignot et al., 2011), the possible contribution of subglacial melting to the mass balance of the ice sheet (e.g. Pattyn, 2010), and locating ice-coring sites to reconstruct ancient climate (Wolff et al., 2006;Jouzel and Masson-Delmotte, 2010).
A recent numerical modelling experiment of the icesheet thermodynamics showed that 55 % of the grounded part of the Antarctic ice sheet is at pressure melting point, though this estimate is hampered by insufficient knowledge of geothermal heat flow (Pattyn, 2010). Satellite observations of ice-sheet surface elevation changes strongly suggest movement of subglacial water over timescales of years or less (Smith et al., 2009;Fricker et al., 2007).
More direct approaches to assess the subglacial environment are airborne and ground-based radar remote sensing (Carter, 2007;Popov and Masolov, 2007;Siegert et al., 2005). In addition, several studies have analyzed the reflectivity of radio waves at the ice base. This reflectivity approach has been applied to several areas of Antarctica (e.g. Bentley et al., 1998;Carter et al., 2009;Christianson et al., 2012;Gades et al., 2000;Jacobel, 2010;Langley et al., 2010;Peters et al., 2005;Wright et al., 2012;Zirizzotti et al., 2012). As for methodology, accurate estimation of basal reflectivity requires calibrated radar sounding data for power detection and reliable estimates of radio-wave attenuation. Attenuation uncertainty introduces uncertainty into the bed reflectivity and thus affects diagnosis of the bed conditions. Recently, Matsuoka et al. (2011) argued that in most cases, based on a one-dimensional ice flow model, variations in bed returned power are dominated by variations in englacial attenuation rather than bed reflectivity. He argued that both the accumulation rate and geothermal heat flux anomalies can interfere with the interpretation. Consequently, analytical radar algorithms that have been widely accepted are likely to yield large uncertainties in the temperate-frozen diagnosis.
In the present paper, we analyse radar returned power to characterize the subglacial environments for a total distance of ∼ 3300 km in Dronning Maud Land (DML), East Antarctica ( Fig. 1). We observe that the returned power decreases as ice becomes thicker, but this relationship is not present at depths greater than a critical thickness value that varies according to location. We attribute this finding primarily to the difference in bed reflectivity: higher reflectivity is caused by subglacial water. Using this feature, we attempt to delineate temperate and frozen beds, and discuss the location of the temperate beds in terms of surface elevation, ice thickness, and locations of ice divides or fast-flowing ice.

Instruments
Three ground-based, pulse-modulated radar sounders were used to gather the data used in this study, as listed in Table 1. Two of the radar sounders have a centre frequency of 179 MHz (referred to hereinafter as the 179-1 and 179-2 radars), and the other has 60 MHz (referred to hereinafter as the 60 radar). Use of these radar systems allows us to investigate frequency and pulse-width (vertical resolution) dependences of the bed-returned power. The 179-2 radar and the 60 radar were used previously (Fujita et al., 1999;Matsuoka et al., 2002). Different pulse widths were chosen depending on the field season (1996/1997 or 2007/2008) or on the initial scientific target of the measurements (internal layers or ice thickness). In order to measure the thickness of thick ice, longer pulses (1000 ns or 500 ns) were chosen. When shorter pulses were sufficient to detect ice thickness, shorter pulses (250 ns or 350 ns) were used. Data with more than two settings were used to cross-check ice thicknesses and to diagnose bed conditions.

Factors controlling the bed returned power
We investigate the bed returned power [P bed ] dB in decibels as a function of ice thickness H and lateral location. Here, the brackets indicate values that are expressed in decibels. After the two-way travel of the electromagnetic waves between the radar and the ice/bed interfaces, the bed-returned power [P bed ] dB is expressed as where [S] dB is the sum of the instrumental factors related to, e.g. transmission power, gains by amplifiers, cable loss, antenna gain, antenna area, and refraction loss; [R bed ] dB is the reflectivity at the bed; [G s ] dB is loss due to geometric spreading of the electromagnetic waves; [L] dB is the loss due to englacial dielectric attenuation; and [B] dB is the loss due to birefringence effects. The term [G s ] dB is given by where ε is the dielectric permittivity of ice. The magnitude of [B] dB depends on the radar frequency, the strength of birefringence in the ice, and the orientation of the antenna. Even if this effect appears in the real data, the data will simply yield accidental (random) minima in [B] dB , which will not cause a systematic bias (Fujita et al., 2006;Matsuoka et al., 2009).
The term [R bed ] dB depends on the dielectric contrasts and roughness between ice and its substrate. Because the dielectric permittivity and conductivity of water (e.g. Ray, 1972) are much higher than those of ice and rocks, Fresnel reflectivity of an ice/water interface is larger than an ice/bedrock interface by 10-15 dB (e.g. Peters et al., 2005). The roughness of the ice/bed interfaces also affects [R bed ] dB : interfaces yield larger or smaller values of [R bed ] dB if the reflectors are smoother or rougher, respectively. Since the reflectivity is affected by roughness relative to the radio-wave wavelength, [R bed ] dB can be radar frequency dependent. The transition can occur over a short lateral distance because the phase transition between solid ice and liquid water should be distinct.
The term [L] dB is a function of the temperature of ice, the amount of impurities within the ice, and the propagation path length 2H . The attenuation coefficient α dB m −1 is primarily a function of ice temperature T , so [L] dB can be written as where z is the depth axis (positive downward). At radar frequencies below approximately 600 MHz, [L] dB is practically independent of frequency in the temperature range of the polar ice sheets (Fujita et al., 2000). The temperature field within the ice depends on the boundary conditions and the internal conditions; Matsuoka et al. (2011) (Haran et al., 2005). The red lines show the survey routes for the ground-based radar sounding discussed herein. The traces include the route of the JASE traverse (Holmlund and Fujita, 2009;Fujita et al., 2011) connecting two deep ice coring sites, namely, Dome Fuji and EPICA DML. The dotted black lines are radar survey routes not used for radar diagnosis in this paper. The locations of major sites along the routes are listed in the Supplementary information as Table A. Labelled legs are listed in Table 2. The thin blue traces represent ice divides that separate drainage basins.
This equation implies that we can indirectly determine the lateral variation of [R bed ] dB if the lateral variation of the depth-averaged attenuation is minimal. Here, we assume that attenuation varies minimally and smoothly in space. This assumption is likely valid, since both ice temperature and chemistry presumably vary smoothly in space.
In this study, we search for locations where [R bed ] dB varies by 10-15 dB, which can be associated with the temperate/frozen boundary of the bed, and radar-frequency dependence of the locations and magnitudes of the local variations of the returned power. Local variations of bed-returned power are more likely caused by bed reflectivity than attenuation.

Initial data processing
The procedures for the initial data processing were as follows. First, echoes from the bed are tracked so that ice thicknesses H is determined based on the two-way travel time (TWT) from the surface to the ice-sheet bed and the speed of radio waves within ice. Second, the maximum power associated with the bed echo is defined as [P bed ] dB . We extracted peak power of the time-series of echoes from the bed. Then, the effects of geometric spreading were corrected so that the geometrically corrected bed returned power [P c bed ] dB could be derived. Our instrumental calibration data show that variations in [S] dB are less than 2 dB, so that variations in [P c bed ] dB arise from variations in [R bed ] dB and [L] dB (Eq. 3). The data of both H and [P c bed ] dB were averaged using a moving average over a horizontal distance of ∼ 0.3 km in order to increase the signal-to-noise ratio. We rejected data from the present analysis when [P c bed ] dB is within 3 dB of the noise floor, since our instrumental calibration data show that the returned power uncertainty is larger near the noise floor.

Study area
Ground-based radar sounding data were collected for a total distance of ∼ 3300 km across Dronning Maud Land (Fig. 1). The primary data were collected during the Japanese Swedish IPY 2007/2008 Antarctic Expedition (JASE traverse) (Holmlund and Fujita, 2009;Fujita et al., 2011). The secondary data were collected in the vicinity of Dome Fuji in 1996 and 1997 (Fujita et al., , 1999. The JASE profiles are about 2800 km long in total, and the other profiles are 500 km long in total. The continental DEM (Bamber et al., 2009) was used to calculate surface and bed elevations from the radar-derived ice thickness data along the route. The survey routes include two deep ice coring sites: Dome Fuji and EPICA DML. Ice drilling found temperate beds at these sites (Motoyama, 2007;Murshed et al., 2007). The survey routes also include four subglacial lakes previously detected using radar data (Popov and Masolov, 2007). For coordinates of the selected sites along these profiles, see Table A in the Supplement.
We defined multiple legs of a few hundred kilometres in length, depending on the dynamic conditions of the ice sheet (Figs. 1 and 2; Table 2), and categorised as plateau (P), midstream (M), and coastal (C) legs. P legs have surface elevations higher than ∼ 3400 m, M legs between ∼ 2500 m and ∼ 3400 m, and C legs below ∼ 2500 m. C legs are within a few hundred kilometres of the coast. For many legs, data with more than two settings were available. Such data were used to investigate the effects of different radio frequencies or different pulse widths. For convenience of discussions in this paper, we tentatively classify regions of Antarctica into the following three groups.
(i) Plateau region in the vicinity of the ice divide: There are nine P legs (P1-P9). Along these legs, the ice thickness H ranges between 1780 m and 3460 m. Except for P9, the legs lie within approximately 150 km of either the Dome Fuji summit or the most distinct ice-flow divide of DML ( Fig. 1). P9 is located up to 350 km away from the ice divide. Ice flow speeds are slower than ∼ 2 m yr −1 along P1-P8 legs (Huybrechts et al., 2009;Motoyama et al., 2008;Seddik et al., 2011) and less than 4 m yr −1 along the P9 leg (Motoyama et al., 2008).
(ii) Midstream regions: The midstream M legs are subgrouped into the western side of the traverse near EPICA DML (MW1 and MW2 legs) and the eastern side of the traverse in the Shirase Glacier drainage basin (ME legs). Bed echoes were not continuously recorded along the MW1 leg so it was excluded from further analysis. Leg MW2 is located between EPICA DML and Site-1 adjacent to the inland mountains. Towards Site-1, the ice becomes thinner from ∼ 2900 m to ∼ 400 m ( Fig. 2) and the ice-flow speed decreases to ∼1 m yr −1 (Wesche et al., 2007;Rybak et al., 2007). Along the ME leg, ice thickness varies by ∼ 1000 m with a dominant periodicity of ∼ 20 km and ice flow speeds are relatively high (up to ∼ 18 m yr −1 ) (Motoyama et al., 1995). These conditions are quite different from the P and MW2 legs.
(iii) Coastal area: The C legs are also subdivided to eastern CE and western CW legs. The CW leg is in the vicinity of nunataks, so ice thickness varies greatly (∼ 2500 m to < 100 m), and across Veststraumen ice stream ( Fig. 1) The Cross-sectional plot of the ice sheet. Labelled legs are discussed herein. The abscissa indicates the horizontal distance along the profiles, and the ordinate indicates the elevation above the WGS84 ellipsoid. The vertical exaggeration is approximately 500 times. The locations of major sites are shown as vertical dashed lines. The bed elevation (red, blue, and green) was calculated from the surface elevation data and ice thickness data extracted from the radar sounding data (Fujita et al., 1999(Fujita et al., , 2011, and where this dataset does not show the bed reflection, the data with thick black dashed lines are from Huybrechts et al. (2009). The red, blue, and green traces are predicted temperate bed, frozen bed, and uncertain (or intermediate) bed conditions, respectively, as the diagnosis result of the present work. The accumulation data were obtained using stake measurements for the section between S16 and Dome F (Motoyama et al., 2008), ice and snow radar surveys between Dome F and EPICA DML (Fujita et al., 2011), and firn-core studies between EPICA DML and Wasa (Rotschky et al., 2007). (Näslund et al., 2000). The CE leg crosses a tributary of Shirase Glacier.

Plateau region legs P1, P2, and P3
The relationships between 179-MHz [P c bed ] dB and H along the P1, P2, and P3 legs are shown in Fig. 3a. We hereinafter refer to this type of graph as an H-P plot. The distribution of the data reveals that [P c bed ] dB tends to decrease with increasing H . At H >∼ 2800 m, [P c bed ] dB tends to either increase or become independent of ice thickness. Here, we refer to this large-scale tendency of the data distribution as similar to a hockey-stick. The noise floor in [P c bed ] dB is independent of the two-way travel time, but with geometric correction the noise floor apparently increases with ice thickness in the H-P plot; the hockey-stick shape is observed 5-10 dB above the noise floor. Ice thickness dependence of 60-MHz [P c bed ] dB for the same legs show a similar hockey-stick shape in Fig. 3c.
The spatial distributions of H and [P c bed ] dB with respect to horizontal distance x for two different frequencies are given in Fig. 3b for the 179-MHz data and in Fig. 3d for the 60-MHz data. We hereinafter refer to this type of graph as an X-PH plot. In the X-PH plot, the scales of the left axis for ice thickness and the right axis for the returned power are adjusted using the gradient of the linear least-squares fit in the H-P plot for the handle of the hockey stick. When ice thickness and returned power curves fall one over the other in the X-PH plot, the returned power follows the regional trend that [P c bed ] dB decreases linearly with ice thickness (i.e. the handle of the hockey stick). When these two curves are not close, the returned power is anomalously high for the ice thickness, corresponding to the toe-heel (plate) of the hockey stick. The ice-thickness range used for the linear approximation in the H-P plots is somewhat arbitrary; for these legs, we fit the data shallower than 2600 m and the uncertainty related to this arbitrary value is discussed later.
Scaling the two ordinates in the X-PH plot is equivalent to very rough corrections for englacial attenuation and allows visualisation of possible anomalies of bed reflectivity and their location along the radar transects, highlighting x locations and H values where subglacial conditions alter.
Linear approximation in the H-P plot implicitly assumes that the attenuation rates vary minimally and have no depth dependence (Matsuoka et al., 2011). In reality, the validity of this assumption largely depends on how the data ensembles are made. We also used data ensembles for a fixed distance and for certain glaciological conditions (e.g. ice flow speed). However, none of these show a significant improvement from the method presented above. In reality, the H-P plot for a single leg can include multiple hockey sticks distributions associated with the attenuation variations within a leg. Nevertheless, Fig. 3a-d clearly demonstrates the bulk feature of the hockey sticks: large anomalies in the bed returned power are found in thicker ice. Data locality represented by lines between data points in the H-P plot show distinct features: data in the vicinity of the stick handle show that, as ice thickness varies locally, the returned power varies as the regional trend (handle) predicts. In contrast, for thicker ice, data in the vicinity of the stick heel and toe show that either returned power varies greatly even if ice thickness varies little along the radar transect or that the retuned power does not vary locally even if ice thickness varies significantly along the radar transect.
The anomalous bed returned power is also shown in a more conventional way, δ[P c bed ] dB versus location x, where δ[P c bed ] dB is defined as the difference between radarobserved [P c bed ] dB and the predicted value using the linear trend in thinner ice (Fig. 3e). We refer to this type of plot as an X-δP plot. Along P1-3 legs, δ[P c bed ] dB values are continuously 10 dB over a distance of more than 150 km in the range x = 850-1000 km, and over more than 30 km from x = 1120 km. δ[P c bed ] dB often have anomalous increases (transitions) of about 10-15 dB at x locations near the increasing point of the H-P plot.
Anomalously high values of δ[P c bed ] dB show similar regional patterns at both 60 MHz and 179 MHz ( Fig. 3b and d). Figure 4 shows relationships between [P c bed ] dB measured at different radar frequencies and different radar pulse widths. These relationships are nearly linear and the data are scattered from the linear trend by a few dB at most. This result shows that the hockey-stick feature is fairly common at two different radio frequencies and pulse widths.

Plateau legs P4-P9
The other inland legs P4, P6, and P7 show similar properties to legs P1 to P3 (Fig. 5a-c). The data show a large-scale hockey-stick-like distribution in the H-P plot and anomalously high δ[P c bed ] dB values occur over a distance of more than 90 km in the range x = 240-330 km, and over more than 150 km from x = 370 km. However, there are much fewer data points that constitute the handle of the hockey sticks compared to the P1-P3 legs (Fig. 3).
Legs P5 and P8, which exactly follow the ice-flow divide, show distinctly different features from the other P legs (Fig. 5d-f). At a given ice thickness, the [P c bed ] dB values vary more than the other P legs and their mean values are roughly 10-15 dB larger than the other P legs. In addition, no apparent hockey-stick feature is identified and the associated leastsquares gradient is poorly defined. To express the X-PH plot for legs P5 and P8, we use the least-squares gradient of the neighbouring legs P4, P6, and P7, which are within 130 km of the P5 and P8 legs. We will discuss the meaning of the data features in the discussion section. least-squares gradient for H < 2400 m. The H-P plot shows that, at ice thicknesses greater than 2400 m, local variations of the returned power are significant although the ice thickness does not vary in the vicinity. This feature corresponds to remarkably large δ[P c bed ] dB values over a distance of more than 90 km in the range x = 660-810 km. [P c bed ] dB often fluctuates with an amplitude of 10 dB and within a few kilometres or less. These features are similar to those observed for the neighbouring P1-P3 legs (Fig. 3).

Midstream regions
The H-P plot along the MW2 leg has a well-defined hockeystick shape. Anomalously high [P c bed ] dB values constituting the toe to heel of the stick appear at ∼ 2500 m (Fig. 6a). In shallower ice, the depth dependencies of the returned power is different for ice thickness ranges of <∼ 1500 m and 1500 to ∼ 2500 m. Since we aim to examine anomalously high bed returned power at greater depths, we fitted the data only for the thicker ice ranging between ∼ 1500 and of anomalous returned power if other ice thickness ranges are used for the fitting; however, alternate fittings give different magnitudes of the anomalous power. Figure 6b shows that the [P c bed ] dB and H profiles lie one over the other. The δ[P c bed ] dB is anomalously higher by 10-15 dB at x locations where H >∼ 2500 m, and above a few subglacial mountains (Fig. 6c). An anomalously high δ[P c bed ] dB value is found at the EPICA DML ice core site, where ice drilling found a temperate bed (Murshed et al., 2007).
Data along the leg ME show different properties (Fig. 6df). For this leg, the hockey-stick shape is ill defined, and [P c bed ] dB does not correlate well with increasing H . The trajectory of the data points in Fig. 6d shows that the bed returned power varies 10-15 dB locally even if the ice thickness remains nearly constant, and that the bed returned power varies little locally even if ice thickness varies significantly. Such features were typically found not in the handle but in the toe and heel of the hockey stick for the P legs. As done for the legs P5 and P8, we used the least-squares gradient obtained along the immediate neighbouring P9 leg to develop the X-PH plot.

Coastal area
The H-P plot for the CE leg shows a well-defined hockey stick shape (Fig. 7a). The linear least-squares approximation gives X-PH and X-δP plots showing increase of δ[P c bed ] dB in inland (Fig. 7b). In leg CE, remarkable increases in [P c bed ] dB occur at some locations (70 km < x < 160 km and x > 190 km), where the ice flow velocities are faster (see the flow velocities in Fig. 7b and Supplement Sect. B). Geographically, this area is located directly upstream of the Soya Coast, where Hirakawa (1997, 2002) observed a landform of bedrock that has been eroded by subglacial melt water some time in the past. For x >∼ 250 km, we observed no features indicating synchronization between the variations of [P c bed ] dB and H . This area is located upstream of the Shirase Glacier (Fig. 1), through which most of the ice in this drainage basin flows to the sea. Again, in the lower right-hand side of the H-P plot, tracks appear as either nearvertical or near-horizontal fluctuations.
The CW leg shows anomalously high bed returned power for ice thicknesses of 1300-1600 m, but not a single hockey stick can be defined (Fig. 7d). We interpret that this feature is made of multiple hockey sticks that have heels at  different ice thicknesses. In leg CW, remarkable increases in [P c bed ] dB occurred along x locations across the Veststraumen ice stream (x = 2625-2725 km) and another ice stream at x =∼ 2560 km. In the lower right-hand domain of the H-P plot, tracks appear as either near-vertical or near-horizontal fluctuations.

General data properties
This analysis has identified several features in the data common from the inland to the coast. In most H-P plots, a representative hockey-stick-like distribution of [P c bed ] dB is defined. The heel of the stick is located at an ice thickness ranging between 2600 and 2800 m in the inland region (P legs), 2400-2500 m in the midstream areas (leg MW2), and 1000-1500 m in the coastal area (legs CW and CE). However, no exact single H value represents each of the partitioned legs. [P c bed ] dB often have anomalous increases (transitions) of about 10-15 dB at x locations near the increasing point. These features are observed regardless of the radar frequency or radar pulse width. Local variability of [P c bed ] dB show two distinct phases. One is that [P c bed ] dB tends to be the value predicted with the linear approximation, so the attenuation rate is regionally uniform. The other is that the data within a given vicinity do not follow this linear prediction: [P c bed ] dB varies greatly even if the ice thickness varies little in the vicinity or vice versa. The latter phase is apparent in the vertical and horizontal data trajectories in the H-P plots. The amplitude of the near-vertical fluctuations are often in the range 10-20 dB.

Primary cause of the anomalously high returned power under thicker ice
Our analysis found anomalously large [P c bed ] dB for thicker ice. Such anomalous features can be caused by anomalous bed reflectivity, englacial attenuation, or both (Eq. 3). Here, we model englacial attenuation in the inland region where horizontal advection of heat is minimal so that a onedimensional heat-flow model can be used to realistically replicate the depth profiles of ice temperature. These can then be used to model englacial attenuation. We do not attempt to estimate the attenuation rates accurately, since the boundary conditions, evolution of the ice sheet, and englacial chemistry are not well known. Rather, we examine its qualitative behaviour.
Following Matsuoka et al. (2011), we modelled depth profiles of ice temperature and then the attenuation coefficient α for ice thicknesses ranging between 1700 and 3500 m (Table 2), surface accumulation ranging between 20 mm yr −1 (ice equivalent) and 60 mm yr −1 (Table 2), geothermal fluxes of 40, 50, and 60 mW m −2 , and a surface temperature of −50 • C. Typical values of geothermal fluxes and a surface temperature in the P legs are taken from Fox-Maule (2005), Shapiro and Ritzwoller (2004), and Comiso (2000).
Integrating the attenuation coefficient α gives the two-way attenuation [L] dB (2), which is shown in Fig. 8 as a function of ice thickness, and is similar to a H-P plot but not accounting for the bed reflectivity. Since attenuation rates are ice-thickness dependent, two-way attenuation does not vary linearly with ice thickness even if the boundary conditions remain the same. Smaller ice thickness dependencies correspond to a temperate bed. In the data analysis above, we linearly approximated the ice thickness dependence of [P c bed ] dB for ice thicknesses ranging between 1700 and 2400 m; an ensemble-mean depth-averaged attenuation rate for this ice thickness range is 4.5 dB km −1 . Figure 8 also shows power loss estimated with this ensemble-mean attenuation rate and with the linear least-squares approximation of all ensembles of the returned power for ice thinner than 2600 m. Figure 8 demonstrates that these two linear approximations yield non-zero anomalous δ[P c bed ] dB solely from the variations, suggesting that δ[P c bed ] dB cannot be immediately interpreted as a proxy of bed reflectivity. However, the heeltoe features of the hockey stick cannot be replicated if the boundary conditions remain unchanged. A considered possibility is that boundary conditions alter in such a way that geothermal flux is low (40 mW m −2 ) for thin ice but high for thick ice (60 mW m −2 ) along a 500 km-long radar transect in the vicinity of the flow divide so that the heel-toe features can be replicated solely with the attenuation variations due to geothermal flux, which we think unlikely. A more likely scenario is that the boundary conditions remain similar but the ice thickness varies. In this case, the heel-toe features cannot be predicted by accounting for the two-way attenuation. Therefore, we argue that the hockey-stick shape is caused by alternation of the bed conditions from frozen underlain thinner ice to temperate underlain thicker ice. However, the possible range of δ[P c bed ] dB that can be caused by attenuation variations largely depends on the variations of geothermal flux and surface accumulation rates, so it is impossible to accurately locate the boundary between frozen and temperate beds.
The bed reflectivity can change with the bed roughness in the scale of radio-wave wavelengths and the substance at the base. Although the effect of the target (bed) roughness is radar-frequency dependent (e.g. Fung, 1994;Ulaby et al., 1986), our data show that the primary properties of the data are independent of radar frequency (Fig. 4a). In terms of wavelengths in ice, we used both 0.94 m and 2.8 m, being different by a factor of 3. This difference in wavelength gave no detectable effects. Based on this fact, we exclude bed roughness as a possible major cause of the larger bed returned power at great depths. . The gray dotted line shows the power loss estimated with the ensemble-mean attenuation rate derived using model estimates for ice thicknesses less than 2400 m. The black dotted line shows the power loss estimated with the attenuation rate derived with the linear least-square method applied to the modelled results for an ice thickness less than 2400 m.

Empirical methods to diagnose bed conditions
Our interpretations of the general data properties described in Sect. 3.5 are as follows. In the H-P plots, the representative hockey-stick-like distribution of [P c bed ] dB is due to superposition of many hockey-stick-like distributions with different threshold H and [P c bed ] dB values for onset of the anomalous increases. The variability of the heel position of the stick is explicable naturally by the variable conditions of the ice sheet. According to the analysis in Sect. 4.1, the gradients of the regression lines do not indicate the attenuation rate within ice. Thus, the regression lines in the H-P plots provide only tentative baselines to highlight the anomalous increases in [P c bed ] dB . The presence of water at the bed is the most plausible explanation for the anomalous increases (transitions) of [P c bed ] dB by about 10-15 dB at the heel positions of the sticks (H-P plot) and at x locations near the increasing point (X-PH plot and X-δP plot). Thus, to delineate temperate/frozen bed conditions, we need to identify each point indicating an anomalous increase of [P c bed ] dB in the real data. In addition, we interpret that the two distinct phases of the local variability of [P c bed ] dB are indicators of temperate bed and frozen bed conditions. Based on the interpretations above, we propose analytical procedures for delineating bed conditions. The overview of the approach is as follows. More details are given in the Supplement Sect. D.
(i) Step 1: the H-P plot diagnosis. The H-P plot diagnosis must be performed for each partitioned region with nearly uniform surface mass balance and ice-flow speeds. We first check if we can find a large-scale hockey-stick-like distribution in the H-P plot. Second, we observe in which of the two distinct phases the local variability of [P c bed ] dB can be classified. (ii) Step 2: the X-PH plot diagnosis and the X-δP plot diagnosis. Anomalous increases of δ[P c bed ] dB by more than ∼5 dB are taken as an indicator of the existence of liquid water at the bed. As a tentative indicator, if the anomalous increases of [P c bed ] dB are more than 5 dB, we diagnose that the bed is likely temperate. If the anomalous increases are 0-5 dB, we tentatively diagnose that the bed is in an uncertain condition or in an intermediate condition. It must be noted that δ[P c bed ] dB values cannot be used as a strongly reliable threshold to delineate frozen or temperate beds because of variability in the attenuation term. The anomalous increases (transitions) of [P c bed ] dB in the X-PH plot and the X-δP plot are checked again in the H-P plot as the heel positions of the sticks. This procedure can locate the most likely boundary between the frozen and temperate beds.
(iii) Step 3: cross-check and collect more circumstantial evidence. If hockey-stick-like distribution is not identified in Step 1, comparison of the data with neighbouring legs may help.
(iv) If the procedures above do not provide convincing diagnosis, we diagnose that the bed is in an uncertain condition or in an intermediate condition.
Practically, our procedure delineates temperate/frozen conditions without modelling attenuation of radio waves within ice. Partitioning of the data analysis by leg is equivalent to grouping data that may have similar boundary conditions. We diagnose bed conditions in each leg separately.

Delineation of temperate and frozen beds
We apply the empirical method discussed above to delineate the location of temperate and frozen beds. Since data properties are different for each leg, we take slightly different steps for each; the diagnostic steps taken for each leg are summarized in Table 3. The results of the diagnosis are shown in several figures: in an X-δP plot for each leg (Figs. 3, 5, 6,  and 7), in H-P plots (Fig. 9), and on the cross-sectional map of the ice sheet (Fig. 2). The assessment of each leg is described as follows.

Plateau region in the vicinity of the ice divide
The legs P1-P4, P6, and P7 are diagnosed as a mixture of temperate and frozen beds, depending primarily on the values of H . The bed is temperate in wide zones where   Fig. 9. The graphs of the H-P plot are expressed as trajectories of the data from one location to another. The predicted bed condition, i.e. temperate bed, frozen bed, or uncertain condition, is given. The ID of the leg(s) and the radar settings are indicated on each panel. The red, blue, and green traces indicate temperate bed, frozen bed, and uncertain conditions, respectively. See text for details.
H >∼ 2800 m. Frozen bed conditions only occur where H <∼ 2800 m in this area. At some locations above subglacial mountains, the bed is diagnosed as temperate even where H is as small as ∼ 2400 m. The presence of this tem-perate bed is qualitatively explicable by the lower surface accumulation rate above subglacial mountains in this area (Fujita et al., 2011) and special ice flow conditions near the ice divide (see below in this section). These conditions decrease advection of cold ice from the surface. As a result, the bed can be warmer than surrounding areas.
Along the exact ice divide, the legs P5 and P8 were diagnosed only by comparing the [P c bed ] dB data with the other P legs. Judging from the [P c bed ] dB values that were generally larger than the frozen beds in the neighbouring legs by up to 15 dB, in addition to the features of the trajectories in the H-P plots, we diagnosed almost all points as temperate.
This generally temperate bed along the exact ice divide is explicable by a special ice flow condition known as the Raymond effect (Hvidberg, 1996;Martin et al., 2009;Raymond, 1983). In the ice-divide zone, the flow solution changes compared to the flank flows. The vertical velocity is smaller than away from the divide, which makes advection of cold ice mass from the surface smaller. Consequently, the thermalinsulation effect of the ice sheet occurs more strongly in the ice divide zone than in flank areas. Therefore, the bed can easily reach the pressure melting point. Note that the Raymond effect makes the bed temperate at melting point but will not strongly affect the temperature gradient within ice. Geothermal heat flux will be consumed by the heat of fusion, and the Raymond effect should appear primarily in the melting rate at the base. Martin et al. (2009) demonstrated that the strength and spatial distribution of the Raymond effect depend on the strength of flank flows, the migration history of the ice divide, basal sliding, and so on. Further examination of the Raymond effect is provided in the Supplement Sect. E.
Likewise, for leg P9 the bed was diagnosed as temperate at locations where H >∼ 2500 m.

Midstream regions
Along leg MW2, beds in the subglacial trenches including the EPICA DML core site are diagnosed as temperate. Temperate bed are found at variable H ranges of ∼ 1500 to ∼ 2400 m. This situation is explained by the transitional dynamic conditions from the ice divide near EPICA DML to the area influenced by the damming effect of the subglacial mountains.
The conditions along leg ME are as follows. In neighbouring legs P9 and CE, locations with H >∼ 2400 m and >∼ 1500 m, respectively, are temperate. Therefore, in leg ME, locations with H >∼ 2400 m are most likely temperate. The majority of the data in Fig. 6d show features often observed for a likely-temperate bed condition. For some locations where H <∼ 2000 m, we cannot find clear features to diagnose the bed condition, and thus classify it as uncertain.

Coastal regions
In leg CE, we diagnose that the bed at x > 100 km is temperate. We hypothesize that subglacial water flow is active near Mizuho, where active ice flow and thinning of the icesheet has been reported (Mae and Naruse, 1978;Motoyama et al., 1995). Considering the ice sheet topography and the bed topography, as shown in Fig. 10, the subglacial water can flow toward the Shirase Glacier mainly through a trough zone between Mizuho and MD170. In this trough zone, the bed was diagnosed as entirely temperate based on features of the tracks (Fig. 7a) in which the data have numerous nearvertical fluctuations in the toe-heel part of the hockey-sticklike distribution.
In leg CW, we found temperate beds at some trough locations underneath active ice flows, as indicated in Fig. 3e and f. However, at some neighbouring troughs at x = 2590 km and 2740 km, the bed was diagnosed as frozen. This contrasting situation implies that H is not a factor determining the bed condition in this leg. We speculate that water supply from upstream determines the presence of water at the bed.

Agreement with known information of subglacial conditions
In Fig. 10, we compare the results of the diagnosis with known information related to subglacial conditions, such as reported subglacial lakes (Popov and Masolov, 2007) and boreholes at DF and EPICA DML known to have water at their base (Motoyama, 2007;Murshed et al., 2007). The reported locations of subglacial lakes (Lake ID 89, 91, 93, and 97) and two boreholes all agreed with the results of our diagnosis as temperate. See Table A in the Supplement for exact locations.

Spatial distribution of the subglacial environments
The subglacial conditions were examined in the plot of H (ordinate) versus surface elevation (abscissa) for the legs in the Shirase Glacier drainage basin, including the Dome Fuji area, and for the legs between Wasa and Dome Fuji, respectively, in Fig. 11a and b. We confirm that relatively low surface elevation and relatively thick ice are favourable for a temperate bed. Qualitatively, this is explainable by the heat balance between advection of a cold ice mass from the ice sheet surface and the geothermal heat flux from the bed. A thicker ice sheet acts as a more efficient thermal insulator. In DML, the bed inland of the ice sheet tends to be temperate. Frozen areas of subglacial high mountains are an exception. In contrast, the bed of the coastal areas tends to be frozen. Areas of fast-flowing ice in subglacial lowlands or troughs are an exception. Thus, we propose a general picture that subglacial water is dominantly produced at the bed of the wide inland plateau and that the water is discharged to the sea primarily through the bed of the subglacial lowlands or troughs. Fig. 11a and b demonstrate that such conditions are common in two contrasting dynamic conditions on both the western and the eastern side of DML. The estimated fractions of the temperate and frozen bed conditions along the legs are listed in Table 3 Fig. 10. Predicted bed conditions are shown on a bed topography map of DML for the same area shown in Fig. 1. The red, blue, and green dots indicate sites of temperate, frozen, and uncertain bed conditions, respectively. The bed topography map is from the BEDMAP compilation (Lythe et al., 2001), which does not reflect the new data. The bed elevation is shown by the colour scale. Surface elevation is shown by thin black contours. Known information related to subglacial conditions is given in the figure.   Fig. 11. Suggested bed conditions are plotted on graphs of ice thickness H (ordinate) versus surface elevation (abscissa). Panels (a) and (b) are for legs between S16 and the DF area and for legs from Wasa to the DF area, respectively. Frozen dry bed and temperate wet bed conditions are indicated as blue and red traces, respectively. Green traces indicate uncertain bed conditions. The thin black lines represent published data (Huybrechts et al., 2009). The locations of two ice coring sites, i.e. Dome Fuji and EPICA DML, are indicated by yellow symbols.
bed, frozen beds and beds with uncertain condition comprise 56 %, 14 % and 14 % of the total leg distance, respectively, not including leg MW1. For the remaining 16 %, ice thickness was undetected with our radar systems. In such cases, we speculate two conditions: either the ice was too thick for the thickness to be detected or bed reflectivity was small due to a frozen bed. The fraction of the temperate bed to the frozen bed estimated here is similar to that predicted by the ice-flow model (55 % temperate, 45 % cold over the whole of Antarctica) (Pattyn, 2010). We provide a direct comparison between the modelling and our diagnosis in the Supplement Sect. F, which provides good qualitative agreement. We argue that the present analysis/diagnosis provides useful boundary conditions of temperature to ice-flow models.
Previous studies have shown that the major ice stream systems have complex systems of tributaries that extend far inland (e.g. Bamber et al., 2006Bamber et al., , 2000Bennett, 2003;Pattyn et al., 2005;Rignot et al., 2011). These studies suggested that as much as 90 % of the discharge from the Antarctic ice sheet drains through a small number of fastmoving ice streams and outlet glaciers fed by catchment areas. Rignot et al. (2011) presented a digital mosaic of ice motion over Antarctica assembled from satellite interferometric synthetic-aperture radar data. They emphasized the importance of basal-slip-dominated tributary flow over deformation-dominated ice sheet flow. In the Supplement Sect. B, we compare the results of our radar diagnosis with the digital mosaic of ice motion in Antarctica compiled by Rignot et al. (2011). We observe that relatively high/low velocities are often correlated with diagnosed temperate/frozen bed, respectively, which supports the argument made by Rignot et al. (2011). Bell et al. (2011) found that, in the vicinity of Dome A, the base of the ice sheet has refrozen ice over 24 % of the surveyed area. The subglacial melting likely occurs at an ice thickness of ∼ 2900 m (Figs. 2 and 3 in Bell et al. (2011)). Their estimate is comparable with our analysis along the P1, P2, and P3 legs. In the vicinity of Dome Fuji, 900 km away from Dome A region, the bed is temperate where the ice thickness is more than ∼ 2800 m. We found no features that can be immediately interpreted as the evidence of the subglacial refreezing-like plume described in Bell's paper. This is probably because large variations in the ice thickness in the Dome A area provide thermal conditions preferred for basal melting and freezing over a short distance, whereas the Dome Fuji area has much smoother ice thickness variations, so the basal melt water drains to downstream areas rather than being kept locally as refreezing ice. Discharge of subglacial water from the inland to the coast is found in the Aurora basin (Wright et al., 2012), and we speculate that the subglacial setting in the Dome Fuji/Shirase basin is more similar to the Aurora basin than Dome A. The threshold ice thickness between temperate and frozen beds is larger at Dome A than Dome Fuji. This regional difference may be caused by geothermal heat flux underneath the Dome A area (Fox-Maule, 2005;Shapiro and Ritzwoller, 2004).

Conclusions
In this study, we have developed an empirical method to delineate temperate and frozen beds in Antarctica. This method examines relationships between ice thickness and bed returned power, spatial patterns of the bed-returned power anomalously higher than the predicted value using a linear approximation of the ice thickness dependence of the returned power, and also fluctuation behaviours of the data. We applied this method to legs roughly 100-400 km in length with similar glaciological conditions such as ice flow speeds and surface mass balance. Despite its simplicity, we argue that this new method presents plausible results. Our analysis shows that at least 56 % and 14 % of the total leg distance have a temperate and a frozen bed, respectively. The beds with uncertain condition comprise 14 % of the total leg distance (Table 3).
Localized water flow in subglacial troughs and lowlands appears to occur due to water production and supply from upstream. Except for such water discharge paths, the remaining areas near the coast are mainly frozen to bedrock.
Future research efforts should include the compilation of a subglacial temperate/frozen condition map, similar to the BEDMAP bed elevation map (Lythe et al., 2001). Ice sheet modelling using the results of the radar diagnosis as the boundary conditions will be useful for improving the understanding of geothermal heat flux. We note that temperate conditions at the bed are favourable for ice thickness detection by radars because the bed reflectivity increases by up to approximately 20 dB. This advantage needs to be recognised in planning future radar surveys of Antarctica.