Alaska's arctic glaciers have retreated and thinned during recent decades, and glaciers in the central Brooks Range are no exception. Digital elevation models (DEMs) reconstructed from topographic maps (from 1970 and 1973) were differenced from a 2001 interferometric synthetic aperture radar DEM to calculate the volume and mass changes of 107 glaciers covering 42 km2 (1970/1973) in the central Brooks Range, Alaska, U.S.A. For each glacier the 1970/1973 DEM was 3-D co-registered (horizontal and vertical) to maximize agreement between the non-glacierized terrains of both DEMs. Over the period 1970–2001, total ice volume loss was 0.69 ± 0.06 km3 corresponding to a mean (area-weighted) specific mass balance rate of -0.54 ± 0.05 m w.e. a-1 (± uncertainty). The arithmetic mean of all glaciers' specific mass balance rates was -0.47 ± 0.27 m w.e. a-1 (± 1 std. dev.). A value of -0.52 ± 0.36 m w.e. a-1 (± 1 std. dev.) was found when 3-D coregistration is performed over the entire domain instead of individually for each glacier, indicating the importance of proper co-registration. Glacier area, perimeter, boundary compactness, mean elevation, and mean slope were correlated with specific balance rates, suggesting that large, low-elevation, elongated and shallow sloped glaciers had more negative balance rates than small, high-elevation, circular, and steep glaciers. A subsample of 36 glaciers showed a mean area reduction of 26 ± 16% (±1 std. dev.) over ∼35 years.
Climate change is impacting Alaska's glaciers, resulting in accelerated rates of mass loss (Arendt et al., 2002; Molnia, 2007; Berthier et al., 2010). In contrast to the large ice masses in Alaska's south and southeast, Alaska's arctic glaciers, geographically defined as located north of the Arctic Circle (66°33′N), include relatively small valley and cirque glaciers. Most of these glaciers started to retreat in the 1890s with a more significant retreat and thinning during the last four decades (Rabus and Echelmeyer, 1998; Nolan et al., 2005; Sikorski et al., 2009). Although small contributors to rising sea level, these glaciers are important indicators of climate change and provide information on long-term climate variations in an area that has few meteorological stations. In addition, these thinning glaciers contribute additional water to streamflow. Recent observations in the Arctic have shown an increase in discharge levels (Peterson et al., 2002). Continued glacier retreat and thinning will ultimately lead to reduced streamflow (Hock and Jansson, 2005), potentially impacting arctic stream ecology (Nolan et al., 2011). Currently the only long-term glacier mass balance record in Arctic Alaska comes from McCall Glacier (69°18′N, 143°48′W), with an area-averaged rate of -0.35 ± 0.07 m w.e. a-1 between 1956 and 1993 and an increased rate of -0.47 ± 0.03 m w.e. a-1 from 1993 to 2002 (Nolan et al., 2005).
Many studies in Alaska have assessed geodetic glacier mass balance by differencing DEMs of different years (e.g. Cox and March, 2004; Larsen et al., 2007; Berthier et al., 2010). Cox and March (2004) compared mass balance results obtained from both direct and geodetic methods. Larsen et al. (2007) found substantial glacier thinning on glaciers in southeast Alaska and Canada using DEMs derived from U.S. Geological Survey (USGS) topographic maps and Shuttle Radar Topography Mission (SRTM). Berthier et al. (2010) estimated that Alaska's glaciers contributed 0.12 ± 0.02 mm to sea level rise between 1962 and 2009, which corresponds to 7.5% of the total sea level rise estimate for the period 1961–2003 by Domingues et al. (2008).
To determine surface elevation change with DEM differencing, accurate co-registration between DEMs is required. This typically includes a 2-D (horizontal) or 3-D (both horizontal and vertical) co-registration between two DEMs. Several studies have co-registered DEMs, but few have compared the impacts on results from the mis-registration of DEMs prior to differencing. Van Niel et al. (2008) found that a mis-registration of half a pixel dramatically impacted elevation differences, with errors compounded on steep slopes. Nuth and Kääb (2011) described the analytical solution to DEM mis-registration. Berthier et al. (2004, 2006, 2010) co-registered DEMs by minimizing the elevation error of the non-glacierized regions. In this study, we co-register DEMs for each glacier individually and collectively within one spatial domain, using statistical minimization methods for comparison.
This study determines the volume and mass change between 1970 and 2001 for 107 small valley and cirque glaciers in the central Brooks Range, Arctic Alaska, using DEMs derived from USGS topographic maps and airborne interferometric synthetic aperture radar (referred hereafter as the 2001 DEM). For each glacier, a 1970 or 1973 DEM was created from USGS topographic maps and 3-D co-registered to the 2001 DEM to maximize agreement between non-glacierized terrains of both DEMs. To explore the effect of mis-registration, we computed glacier volume change for alternate co-registration methods. We also determined volume change using the 1970/1973 USGS national elevation data set (NED) for comparison to the results derived from our reconstructed DEMs. Optical satellite remote sensing (Quickbird; Digital Globe, http://www.digitalglobe.com) was used to update glacier boundaries and determine area change rates between 1970/1973 and 2001–2007 for a subset of 36 glaciers visible in imagery. Additionally, several topographic and geometric indices were related to mass balance results to explore the dependence of small glacier response on topography to changing climatic conditions.
Alaska's arctic glaciers are located in the Brooks Range, a 1000-km-long mountain range in northern Alaska with elevations reaching 3000 m above sea level (a.s.l.) in the east. The glaciers evaluated in this study fall within the Endicott Mountains subrange (Fig. 1). The glaciers are located at a mean elevation of 1747 m a.s.l. (1250–2210 m a.s.l.). The 107 glaciers investigated here (total area = 42 km2, 1970/1973) represent 13% of ∼850 Brooks Range glaciers covering 520 km2 (the value is based on digitized glacier outlines from USGS maps; B. Manley, unpublished data) and corrects a previously published erroneous area estimate (598 km2) by Berthier et al. (2010) (Berthier, personal communication, 2011). The majority of the 107 glaciers are small, north-facing cirque glaciers (Fig. 2) ranging in size from 0.05 km2 to 1.97 km2 with mean size of 0.39 ± 0.35 km2 (±1 std. dev.). The mean annual temperature and annual precipitation totals for the period 2000–2009 observed at the Atigun Pass weather station were -9.5 °C and 0.58 m, respectively (Snowtel Site #957, 1463 m a.s.l., NWCC http://www.wcc.nrcs.usda.gov/snotel/Alaska/alaska.html; Fig. 1).
The most recent DEM (10 m grid spacing) was acquired 18–30 August 2001 by the Intermap Technologies Corporation's Star3i system. Two separate areas covering over 19,000 km2 were captured and contained 123 glaciers within the boundaries (Fig. 1). A total of 16 glaciers were eliminated from the analysis—nine due to small size (<0.05 km2) and seven due to poor agreement in non-glacierized elevations between the earlier and the 2001 DEM. The proprietary Star3i system uses interferometric synthetic aperture radar (IFSAR, X-band) to generate DEMs with a horizontal accuracy of 2 m and vertical accuracy of 1 m (Intermap, http://www.intermap.com). Wave penetration into snow from X-band is <3 cm when water is present (Haritashya and Singh, 2011). The 2001 DEM was collected during the late summer period when ablation areas are generally snow-free and water is present in the snow in the accumulation area, so wave penetration is considered minimal in this study.
NATIONAL ELEVATION DATA SET (NED) 1970/1973 DEM
The 1970/1973 USGS national elevation data set (NED) DEM (∼40 m grid spacing) was downloaded from the USGS web server ( http://seamless.usgs.gov) and projected into UTM zone 5 WGS84 along with re-sampling to a 10 m grid spacing (cubic convolution method). The NED DEM was originally created from vertical photography derived elevation contours produced for USGS topographic maps. The specific creation methods are poorly documented, and contour interpolation techniques are unclear. These unknown aspects of the NED prompted us to construct new DEMs from the original USGS topographic maps with the same grid spacing as the 2001 DEM (10 m).
RECONSTRUCTED 1970/1973 DEMs
We reconstructed the 1970/1973 DEMs using a combination of manual and automated work flows, henceforth referred to as the reconstructed DEMs. Eleven original 1:63,360 USGS topographic map separates (mylar sheets with brown lines for elevation contours and blue lines for glacier boundaries) were acquired digitally (1200 dpi or 1.3 m ground resolution) from the USGS. The maps referred to 1970 air photos for 102 glaciers and 1973 air photos for five glaciers. An example of a glacier's 1970 NED DEM, reconstructed 1970 DEM, and 2001 DEM is illustrated in Figure 3, pails a—c. Within a Geographical Information System (GIS, ArcGIS v9.3), map separates were georeferenced (mean georeferencing error: 2.2 ± 1.0 m) using corner coordinates. Glacier outlines were digitized from the original USGS maps.
Contours within each glacier and within a 1 km distance of each glacier boundary were automatically digitized from georeferenced map separates and manually attributed with elevation values. Contours were interpolated to reconstruct DEMs using a modified spline technique with drainage enforcement (Anudem v5.2) root mean square error (RMSE = 2.4 ± 0.2 m, mean ± 1 std. dev.) at the same 10 m grid spacing as the 2001 DEM (Hutchinson, 1989). Contour interpolation with no drainage enforcement yielded no difference within results.
Due to the small glacier sizes within the study, proper coregistration of DEMs was necessary prior to DEM differencing. A preliminary investigation indicated a need for 3-D co-registration (both horizontal and vertical) of reconstructed DEMs to match the more accurate 2001 DEM, thereby reducing potential inaccuracies in volume change calculations from mis-registration. We assumed that the 2001 DEM represents the true surface elevation and thus treated it as a reference DEM. Its non-glacierized area was assumed to have experienced no surface elevation change over the investigated time period. For each glacier, the 1970/1973 reconstructed DEM was shifted sequentially by up to 10 pixels in all combinations of the 4 cardinal directions. A root mean square error (RMSE) was used to estimate the vertical error between the reconstructed DEMs and the 2001 DEM. RMSE values were only calculated for the non-glacierized terrain surrounding the glaciers (1 km buffer). Each individual reconstructed glacier 1970/1973 DEM was co-registered to a final location based on the minimum RMSE values found for the non-glacierized region. Each DEM was then vertically shifted by ± 25 m at 0.1 m intervals to allow 3-D co-registration. A RMSE value was again calculated for each individual vertical shift combination between each individually reconstructed DEM and the 2001 DEM. The initial 2-D co-registration of DEMs reduced the associated horizontal positional error due to georeferencing errors. As of 2011, no vertical datum transformation exists in Alaska for NGVD29 (USGS maps) to a more commonly used datum (Nolan et al., 2005). The 3-D co-registration eliminated the need for vertical datum transformations from NGVD29 to NAVD88 and reduced concern over the poor survey control of original USGS maps. We investigated the validity of the two-step procedure adopted here for co-registration by performing co-registration in one single step for 20 glaciers. Differences in results were negligible.
For comparison we also computed mass balances without any DEM co-registration and with only a vertical but no horizontal coregistration. We also applied all three co-registration techniques to a single DEM containing all glaciers (individual reconstructed DEMs aggregated into one DEM), i.e. the entire DEM was co-registered rather than each glacier individually. The purpose was to assess whether or not this simpler method yields results similar to those obtained from co-registering each glacier DEM individually. Finally, we applied the same co-registration methods to each individual NED glacier DEM. RMSE and co-registration parameters are shown in Table 1. As expected, the elevation error of the non-glacierized areas decreased as co-registration was performed in more directions. Figure 4 visualizes the direction and magnitude of the 3-D co-registration of each individual glacier DEM and indicates large variability in direction and magnitude between glaciers.
MASS BALANCE CALCULATION AND GLACIER AREA
Volume change was calculated for each glacier by subtracting the reconstructed 1970/1973 DEMs from the 2001 DEM. We assumed lost volume consisted of ice (density = 900kgm-3; Bader, 1954) to allow volume change conversion to mass change. Hence, our estimates are an upper bound since the actual density may be lower due to retreat and thinning of firn; however, information to quantify this effect is lacking. Specific mass balance was determined by dividing the volume change by the glacier area. Previous studies have generally used the mean glacier area over the considered period for calculating specific mass balances (Finsterwalder, 1954; Arendt et al., 2002). For the earlier DEMs we determined glacier area from the glacier outlines digitized from the maps. Our 2001 glacier area estimate is based on manual digitization of glacier boundaries on recent high-resolution optical satellite images (Quickbird) via Google Earth™ (Fig. 3, part e). Of the 107 glaciers, 20 were not visible in satellite images (covered by snow, clouds, or shadows). All other glaciers visible in satellite images showed signs of retreat. A total of 8 out of the 87 remaining glaciers had no ice present, but had past glacier evidence (e.g. defined lateral or terminal moraines). Several glaciers were visible on satellite imagery, but glacier boundaries were difficult to discern (e.g. due to rockfall). We were able to digitize the glacier outlines of 36 glaciers that had clearly defined boundaries on the imagery. These images were captured in 2001 (1 glacier), 2003 (9 glaciers), and 2007 (26 glaciers).
Co-registration parameters for the 107 investigated Brooks Range glacier DEMs using two different 1970/1973 DEMs (‘reconstructed’ and NED DEM) and several co-registration methods. Results are shown for the case where 3-D co-registration was applied to each glacier individually and for the case where it was applied to the entire DEM domain. ‘Reconstructed’ DEMs were produced in this study from the original U.S. Geological Survey topographic maps. RMSE is the root mean square of the elevation error of the 1970/1973 DEM given as the arithmetic mean of all glacier DEMs ( ± 1 std.dev.). Easting/Northing shift reflects the horizontal directional shift. Positive values indicate a North or East shift, while values reflect a South or West shift ( ± 1 std.dev.).
For each of the 36 glaciers a mean area reduction rate between the earlier (1970 or 1973) and the later (2001, 2003, or 2009) date was calculated. These rates were then applied to each of the 36 glaciers to compute the area in 2001 and 1970 (if glacier outline referred to 1973). The arithmetic mean of all 36 area reduction rates was applied to all remaining glaciers to compute the 2001 area of each of these glaciers. The computed total area of all 107 glaciers for the year 2001 was 30.5 km2. The mean of each glacier's area in 1970 and 2001 was then used to compute the specific balances of each glacier. An average mass balance rate (m w.e. a-1) was calculated for each glacier over the 31- (1970–2001) or 28-year period (1973–2001).
To compute an area-weighted specific glacier mass balance rate for the entire spatial domain for 1970–2001, we extrapolated the 1973–2001 volume losses of the 5 glaciers with 1973 DEMs to the earlier date (1970) by assuming their specific mass balance rates over 1970–2001 to be constant. To convert into specific mass balance units we divided the volume change by the mean of the 1970 and 2001 glacier areas.
TOPOGRAPHIC AND GEOMETRIC INDICES
Within a GIS, we used each glacier's boundary from 1970/ 1973 and the 2001 DEM to calculate glacier area, perimeter length, mean elevation, mean slope, boundary compactness, and mean annual potential solar radiation. These were calculated to assess if these factors impact the response of each glacier to changing climate. Boundary compactness is the ratio of the perimeter of a circle with the same area as the glacier and the perimeter of that glacier (Allen, 1998). Boundary compactness ranges between 0 and 1, where elongated polygons typical of valley glaciers have lower values than more circular polygons typical of cirque glaciers. Mean annual potential solar radiation was also calculated for each glacier using topography from the 2001 DEM. We tested indices for significant (p < 0.05) Pearson correlations (SPSS v18) with specific mass balance rates for all glaciers. Because several indices are naturally correlated (e.g. greater area results in greater perimeter), principal component analysis was also conducted to identify a set of uncorrelated variables.
Standard principles of error propagation were used to determine error estimates between two DEMs within the DEM differencing process (Burroughs et al., 1998). However, DEMs are highly spatially autocorrelated. To account for the spatial autocorrelation, we estimated the uncertainty using uncorrelated measurements between the two DEMs (Nuth and Kääb, 2011). We assumed an autocorrelation distance of 0.25 km to allow for an adequate sample size. Other studies' autocorrelation distances range from 0.1 to 1 km (Koblet et al., 2010; Kääb, 2008). Thus, we examined 4% of the total population of pixels differenced between DEMs to account for spatial autocorrelation. We estimated the uncertainty between the two independent DEMs by the following equation:
where σ1 and σ2 represent the random errors associated with each individual DEM. These were derived from the error components listed in Table 2. These include the original glacier contour representation on original USGS maps and contour interpolation. We calculated the standard error following the methods of Nuth and Kääb (2011) resulting in an overall uncertainty of 0.05 m w.e. a-1 for the area-averaged specific balance of all 107 glaciers. This standard error value is low in comparison to the individual glacier errors (ranging from 0.08 to 0.13 m w.e. a-1 mean = 0.11 m w.e. a-1) due to the larger number uncorrelated measurements.
Sources of random error for each individual pixel of the DEMs.
The error associated with the original USGS topographic maps dominates the uncertainty in our reconstructed DEMs. Aðalgeirsdóttir et al. (1998) found errors as great as ± 45 m in glacier accumulation areas of the Harding Icefield for USGS contour elevations. However, errors tend to be larger in wide accumulations areas with small slope (Aðalgeirsdóttir et al., 1998; Muskett et al., 2003). We assumed a smaller error of ± 30 m for the USGS contour elevations as our investigated glaciers are narrow valley or cirque glaciers. The 2001 DEM (vertical accuracy ± 1.0m) was collected using X-band IFSAR during the late summer period. We assumed minimal penetration of radar waves into glacier snow/ice and so allowed a 0.1 m error into our error estimates. Any other errors were considered small compared to the elevation change errors and not considered here.
Mean mass balance rates, , for the 107 investigated Brooks Range glaciers for the period 1970/1973 to 2001 using two different 1970/1973 DEMs and several co-registration methods. ‘Num. +’ is the number of glaciers found to have mass gain along with the percentage of total number of glaciers in parenthesis. is the arithmetic mean of the specific balance rates of all glaciers ( ± 1 std.dev.). min /max is the minimum/maximum mass balance rate of the glacier sample (m w.e. a-1). ΔZ is the mean elevation change of all glaciers (± 1 std.dev.).
Results and Discussion
MASS BALANCE 1970–2001
Using 3-D co-registered reconstructed 1970/1973 glacier DEMs, the total ice volume lost by the 107 glaciers between 1970 and 2001 was 0.69 ± 0.06 km3, corresponding to a mean specific area-weighted balance rate of -0.54 ± 0.05 m w.e. a-1 ( ± uncertainty). The arithmetic mean of all balances (not area-weighted) is - 0.47 ± 0.27 w.e. a-1 (± 1 std. dev.; Table 3, part a). The rate varied widely among the glaciers ranging from -1.12 ± 0.09 (± uncertainty) to 0.34 ± 0.12 m w.e. a-1 with 104 glaciers losing mass and three glaciers gaining mass (0.05 ± 0.12, 0.09 ± 0.12, and 0.34 ± 0.12 m w.e. a-1) during the study period (Table 3, part a). These three glaciers are among the smallest glaciers (0.09, 0.09, and 0.05 km2, respectively). Hence, the positive balances may be artifacts of co-registration errors.
The observed volume loss is consistent with an observed increase in air temperatures and a decrease in annual precipitation. Weather stations are sparse in the region; only two long-term records are available for the Arctic—Barrow and Barter Island—and are roughly 300 and 100 km away from the studied glaciers, respectively (Fig. 1). Mean temperature increased at both Barrow and Barter Island between 1965 and 1995 (1.4 °C and 1.0 °C, respectively; Curtis et al., 1998) while annual precipitation decreased between 1949 and 1988 (30% and 47%, respectively; Curtis et al., 1998). On McCall Glacier, Rabus and Echelmeyer (2002) estimated a 1.1 ± 0.3 °C increase in air temperature between 1975 and 1995.
Our mass losses over the period 1970 to 2001 (0.54 ± 0.05 m w.e. a-1) are larger than those from previous studies within the Brooks Range. McCall Glacier, located ∼220 km from our closest investigated glacier (Fig. 1) showed mass balances of -0.35 ± 0.07 m w.e. a-1 from 1956 to 1993 and -0.47 ± 0.03 m w.e. a-1 from 1993 to 2002 (Nolan et al., 2005). The glacier is much larger (∼6.4 km2) than our glaciers. Berthier et al. (2010) reported a mean mass balance rate of -0.37 ± 0.06 m w.e. a-1 for all Brooks Range glaciers during 1956–2001. They applied McCall Glacier's mass balance profile to the area-altitude distribution of the unmeasured glaciers based on the conclusion of Rabus and Echelmeyer (1998) that McCall Glacier was representative of glaciers in the northeastern Brooks Range. Our results indicate that McCall Glacier's mass changes may not be representative for the regional-scale mass changes occurring in the Brooks Range, although differences may at least partially be due to different time periods between our and their study.
Despite the arctic location of our investigated glaciers, our results are consistent with those from glaciers located throughout Alaska. Arendt et al. (2002) found an area-weighted mean thinning rate of -0.52 m ice eq. a-1 for 67 glaciers located throughout Alaska between the mid-1950s and mid-1990s. Berthier et al. (2010) in their statewide Alaska DEM differencing study found varying regional results for the 1960s–2006 time period ranging from -0.19 m w.e. a-1 in the Alaska Peninsula to -0.65 m w.e. a-1 in southeast Alaska. Although these time periods are different, all studies show consistent thinning over the last 50 years of the same order of magnitude.
GLACIER AREA CHANGES
For the 36 glaciers whose recent boundaries could be digitized, total glacier area was reduced from 21.5 km2 (original USGS maps from 1970 or 1973) to 16.4 km2 (satellite images between 2001 and 2007), a 24% area reduction in about 35 years. The arithmetic mean of all individual glaciers' area reduction rates was 0.73 ± 0.44% a-1 (± 1 std. dev.). Results may be affected by errors due to misclassification of perennial snow in the USGS maps; this error is unknown. However, Granshaw and Fountain (2006) compared a glacier inventory derived manually from aerial photos (Post et al., 1971) to a new digital inventory derived from USGS topographic maps. The total glacier-covered area of the former (116 km2) differed from the latter by only 1.5 % over a total glacier area of 116 km2. These 321 glaciers were very similar in area (mean area = 0.4 km2) to our investigated glaciers (0.39 km2).
Figure 5 depicts annual reduction rates versus original 1970/ 1973 glacier area. Despite the scatter, results indicate that relative area loss tends to be greater for smaller glaciers. Results may not be representative of the entire glacier sample since the size distribution of the 36 glacier subsample is biased towards larger glaciers (Fig. 2). Nevertheless, these results are an additional indication of widespread mass loss and are consistent with our quantitative results obtained from DEM differencing.
We tested for correlations between the specific mass balance rates and topographic indices using the Pearson correlation coefficient. Results are based on 3-D co-registration of each individual reconstructed DEM. Table 4 provides statistical summaries of the topographic indices calculated. We observed significant negative correlations between mass balance and glacier area (r = -0.33; Fig. 6, part a) and glacier perimeter (r = -0.39). This suggests that the larger glaciers thinned more than the smaller glaciers. The larger glaciers in Brooks Range are valley glaciers while the smaller glaciers tend to be more circular cirque glaciers. Significant positive correlations were found for mean glacier slope (r = 0.48; Fig. 6, part c), mean glacier elevation (r = 0.19; Fig. 6, part b), and compactness (r = 0.26; Fig. 6, part d; Table 4). This is consistent with the findings above since the small glaciers tend to be steeper, at higher mean elevations, and more circular than the larger glaciers. No significant correlations were found between mass balance rates and potential solar radiation, longitude, or latitude (Table 4).
As these indices are naturally correlated, we used principal component analysis (varimax rotation method; IBM, 2011) to reduce six variables (size, perimeter, compactness, elevation, slope, and radiation) into two components. The first component axis was principally determined positively by glacier size (loading = 0.89) and perimeter (0.95), while negatively by compactness (0.73) and mean slope (0.65). This “shape and steepness” axis represented a gradient from large, elongate, and flat glaciers to small, round, and steep glaciers and captured 49% of the variance among the six variables. The mass balance rate was significantly and negatively correlated with the “shape and steepness” axis component (r = - 0.45, P< 0.01, n = 107), confirming that larger, more elongated, and shallow-sloped (valley) glaciers thinned more than smaller, more circular, steeper-sloped (cirque) glaciers. The second component, which accounted for 23% of the variance, included mean elevation (0.89) and potential solar radiation (0.78). This “elevation and solar radiation” axis represented a gradient from low-elevation glaciers with low potential solar radiation to higher elevations with high solar radiation. The mass balance rate was not significantly correlated with the “elevation and solar radiation” axis component (r = 0.07). Our results are consistent with findings by De Beer and Sharp (2009), who evaluated glacier area retreat for small glaciers (<0.4 km2) in the Monashee Mountains, British Columbia, Canada. They found that larger glaciers exhibited larger retreat rates than smaller glaciers as the smaller glaciers tended to occupy topographically sheltered sites. The authors suggested that shading and topography reduced mass loss of small cirque glaciers. A dependence of (areal) retreat rate on glacier size has also been found for a larger range of glacier sizes in glaciers in the European Alps (Paul et al., 2004), Norway (Andreassen et al., 2008), and the Canadian Rocky Mountains (Jiskoot et al., 2009).
Pearson correlation matrix between mean specific mass balance rates, (using 3-D co-registered reconstructed DEMs) over the period 1970/1973 to 2001, and topographic and geometric indices. The top three rows are the mean ± std.dev. ( + σ), minimum (Min), and maximum (Max) of the results for the 107 investigated glaciers. A is the glacier area, P is the glacier perimeter, Z is the mean glacier elevation, α is the mean glacier slope, I is mean potential solar radiation, C is boundary compactness, λ is latitude, and is longitude. Bold values indicate significance at the 95% confidence interval.
SENSITIVITY OF MASS BALANCE RESULTS TO CO-REGISTRATION
We also computed mass balance using three co-registration methods of the 1970/1973 reconstructed glacier DEMs. These included (1) no co-registration, (2) no horizontal co-registration and only a vertical shift, and (3) a 3-D co-registration (both horizontal and vertical). The specific mass balance rates (mean of all glaciers) differed by up to 0.12 m w.e. a-1 between the three methods when each glacier was co-registered individually (Table 3, case A, a—c). The standard deviation and range of values tended to increase as co-registration was performed in fewer directions. Hence, although the mean balance rate averaged over all glaciers changed relatively little between different co-registration methods, the rate of individual glaciers was affected considerably as indicated by the larger range in balance rates (Table 3).
Past studies focused on DEM co-registration of larger domains including a large number of glaciers (Larsen et al., 2007; Berthier et al., 2010). Larger DEMs may contain internal distortions introduced from georeferencing or collection methods. We compared the results above with those derived from 3-D co-registration of a single DEM containing all glaciers (Table 3, case C). The elevation error (see RMSE in Table 1, case C) was consistently higher for all co-registration methods compared to the results derived from co-registration of each individual glacier DEM (Table 1, case A). This is not unexpected considering the lack of a consistent regionalscale pattern of co-registration direction and magnitude that was necessary to minimize the elevation error of each individual glacier (Fig. 4). Some co-registration vectors exhibit some level of spatial autocorrelation indicating that for these glaciers co-registration at larger spatial scales could be sufficient. However, it is unclear how to identify these scales a priori.
The mean mass balance rate (averaged over all glaciers) differed by up to 0.27 m w.e a-1 between the two methods (Table 3, cases A and C) indicating that co-registration should be performed for each glacier individually rather than adopting the common method of 3-D co-registering a larger domain including many glaciers. Figure 7 illustrates the scatter between the glacier's specific mass balances using individual 3-D DEM co-registration and those derived from the alternative methods.
NED DEM VERSUS RECONSTRUCTED DEMs
We calculated mass balances by subtracting the 1970/1973 NED DEMs from the corresponding 2001 DEMs to compare to the mass balance results derived with our reconstructed DEMs. Although both DEMs are based on the same data source (i.e. 1970/ 1973 USGS topographic map contours), the mass balance results differed slightly from each other (Table 3, cases A and B). Consistent with the results using the reconstructed DEMs, we found a narrower range in mass balance rates and reduced standard deviations for 2-D co-registration versus no co-registration, and further reduction of range and standard deviation using 3-D co-registration (Table 3, case B). The slight differences may be explained by the unknown original interpolation technique used to create NED versus our modified spline technique (Hutchinson, 1989). Another reason for different results may be the initial difference of grid spacing between NED and our reconstructed DEMs. The NED, originally at ∼40 m grid spacing, was re-sampled to the grid spac- ing of our reconstructed DEM (10 m). With similar elevation errors and mass balance results found between the NED glacier DEMs and our reconstructed glacier DEMs, we conclude that the NED DEM data set is appropriate for use for geodetic mass balance calculations in this region.
DEMs reconstructed from historical USGS topographic maps (1970/1973) differenced from a 2001 DEM within the central Brooks Range, Alaska, allowed us to calculate the volume change and mean specific mass balance rate for 107 glaciers over the period 1970/1973 to 2001. The area-weighted mean balance rate was - 0.54 ± 0.05m w.e. a-1 ( ± uncertainty). Hence, the mass loss is somewhat larger than the one reported for McCall Glacier ( -0.35 ± 0.07 m w.e. a-1 [1956–1993] and -0.47 ± 0.03 m w.e. a-1 [1993–2002]), which is located ∼220 km from our closest investigated glacier (Nolan et al., 2005). The mass loss is consistent with observed increasing air temperatures and decreasing annual precipitation at northern Alaska coastal weather stations during the study period. The correlation of multiple topographical indices with specific mass balance suggests that the larger valley glaciers thinned more than the smaller cirque glaciers that have already retreated into topographically more sheltered sites.
Proper co-registration prior to DEM differencing was found to be important for the 107 studied glaciers that are small, located in steep terrain, with a predominately northern aspect. Co-registration direction and magnitude varied largely between glaciers in the investigated domain. In this study it was necessary to co-register each glacier individually instead of collectively within one larger spatial domain containing all glaciers. This is especially necessary when working with historical USGS-derived DEMs that are known to have poor horizontal and vertical control. Mass balances derived from the reconstructed 1970/1973 DEMs differed only slightly from those obtained using the NED DEMs. This demonstrates that for our study region the existing NED DEMS are appropriate for volume change calculations.
We thank B. Manley (INSTARR, University of Colorado Boulder) for providing a subset of digital glacier boundaries used within analysis. The Star3i DEM was funded by the NASA Commercial Remote Sensing Program. Analysis was supported by grants NSF OPP-0714045 and NASA NNX11AF41G. This work was also supported by the U.S. National Park Service's Arctic Network as part of their Inventory and Monitoring Program (H9910030024/J9910324509 and H9910080028/J9940090220). We also thank A., Arendt, M. Loso, C. Mulder, R. Dial, E. Creely, J. Williams, and classmates D. Podrasky, B.Truessel, and M. Engram for comments that greatly improved manuscript readability.