Study Area

Skagit River drains about 8000 km² of the rugged North Cascades mountains, where the bedrock is highly variable, but is composed mainly of metamorphic and igneous rocks (Tabor and Haugerud 1999). The upper part of the basin from the Cascade River mouth to the east shore of Ross Lake is known as the Crystalline Core because it is composed primarily of gneiss and granite (Figure 2). West of the Cascade River mouth, the basin is dominated by weaker metamorphic rocks schist and phyllite, and volcanic rocks from Mount Baker (3285 m) and Glacier Peak (3213 m). Multiple Pleistocene ice ages brought both extensive alpine and continental styles of glaciation to the Skagit watershed and created the present drainage pattern and topography with 2500 m of local relief.

The watershed is defined by three main regional hydrologic divides that have a strong influence on climate and glacier distribution: the Pacific, North Cascades, and Skagit crests (Figure 2; Riedel et al. 2007). Trending northwest through the center of the watershed is the Skagit crest, which once divided the Skagit into separate Puget Sound and Fraser River bound streams. Continental glaciation breached this divide at Skagit Gorge, and resulted in the capture of the upper Skagit by the westbound lower Skagit, adding significantly to the basin's size (Riedel et al. 2007, Simon-Labric et al. 2014). Skagit crest is breached only at Skagit Gorge, and remains a key hydrologic divide within the watershed, capturing moisture that otherwise would reach the Pacific crest downwind.

Figure 2.

Skagit River watershed on 90m DEM base, showing major hydrologic divides as described by Riedel et al. 2007. Glaciers discussed in the text are underlined. Strreamflow gaging stations on the Skagit River at Concrete (C), and on tributaries Sauk (S) and Baker (B) rivers.

Highland climate in the North Cascades is influenced by both maritime climate on the west and continental climate on the east. Strong gradients of temperature and precipitation occur with distance east from the Pacific Ocean and higher elevation (Mass 2008). The western slopes of the North Cascades are among the wettest temperate mountains on earth, with mean annual precipitation exceeding 6 m on some peaks (Davey et al. 2006). Nearly 75% of the precipitation falls as snow at the altitude of glaciers between October 1 and March 1 (Franklin and Dyrness, 1974). Accumulation ranges from four-to-nine meters of water equivalent (w.e.) on four western glaciers, where it is typically 20–30% higher than it is in the eastern part of Skagit basin (Riedel and Larrabee 2011). Lower snowfall and warmer summer temperatures cause the mean glacier elevation to rise 7 m/km from west to east with increased aridity across the North Cascades (Post et al. 1971).

Air temperature also varies systematically with distance east and elevation. East of Concrete, climate becomes more continental, with ever greater temperature extremes and more of the precipitation falling as snow. Mean annual temperature rises from valley floors at 4.7 °C/km, but lapse rates vary from west (maritime climate) to east (continental), and with aspect and season (Porter 1977, Minder et al. 2010).

In 1959, the largest glaciers in the watershed were found along the Pacific crest south of Cascade Pass, along the Skagit crest north of Cascade Pass, on Mount Baker (17.2 km²), and on Glacier Peak (15.2 km²; Figure 2; Post et al. 1971). Boston Glacier, at the headwaters of Thunder Creek on the Pacific crest, was the largest in the watershed at 7 km², and there were nine other glaciers greater than 3 km². Moisture penetrates deeply into the mountains through Skagit Gorge, as illustrated by the unusually large size of Nohokomeen Glacier (1.8 km²), located well east of the Skagit crest on Jack Mountain. Skagit River has an average daily flow of 0.038 km³ at Mount Vernon (Figure 2; Drost and Lombard 1978). The river's headwaters rise in southern British Columbia, and its largest tributaries enter in the lower half of the watershed and include the Sauk (1896 km²) and Baker (770 km²) rivers.

Data

Data used in this study includes basin-wide measurement of glacial extent in 1959 and 2009, a 10 m DEM, seasonal surface mass balance observations on five glaciers, and basin runoff at three gaging stations operated by the US Geological Survey (USGS; Figure 2). Glacial extent in 1959 was obtained from two sources, including a landmark baseline inventory of glaciers in the North Cascades and digital copies of 1:24,000 scale topographic maps. Post et al. (1971) counted all perennial ice bodies larger than 0.1 km² and measured the area of individual glaciers on 1:38,000 scale topographic base maps. This inventory was aided by the use of oblique aerial photographs and on-site observations. Since the Post et al. (1971) inventory did not develop digital glacier outlines, we used digital copies of 1:24,000 scale topographic maps published by the USGS in 1982. The 10 m DEM was produced as part of the National Elevation Dataset and was based on the same ortho-photographs taken in the late 1950s (herein referred to as 1959).

Individual glacier outlines in 2009 were initially determined by Dick (2013) for the North Cascades region based mainly on 2006 and 2009 National Agricultural Imagery Program (NAIP) ortho-photographs with a horizontal resolution of ± 1 m. We focused on the 2009 NAIP imagery, field visits, and 1:12,000 scale stereo-photographs to re-examine these glacier outlines and create a new geospatial dataset for the Skagit Watershed. The 2009 images were more useful because it was an extremely negative balance year that followed seven consecutive negative years (Riedel and Larrabee 2011), allowing for more accurate identification of glacier margins. Modem outlines for the large glaciers on Mount Baker were compiled by Brown (2011) from the 2006 NAIP imagery, oblique air photographs, and site visits.

Seasonal mass balance observations were made by the USGS and the US National Park Service (NPS) on five glaciers in and near the watershed (Figure 2). Glacial mass balance terms follow Cogley et al. (2011). The methods used at all five glaciers follow a stratigraphic system where measurements were made on a water year basis at times of successive minimum (late September) and maximum (late April) balances (Mayo et al. 1972, Ostrem and Bruggeman 1991). Methods are described further by Riedel et al. (2008) for North Klawatti, Silver, Noisy, and Sandalee glaciers (Figure 2). Krimmel (1996), Bidlake et al. (2010) and Josberger et al. (2007) provide summaries of data collected at South Cascade Glacier. Monitoring by the USGS began in 1957 (Meier 1958), while the NPS monitoring of the other four glaciers has been continuous since 1993 for three glaciers and since 1995 for Sandalee Glacier. Granshaw (2001) extended the mass balance record for these glaciers from 1993 back to 1959 by scaling with the cumulative balance of South Cascade Glacier and geodetic estimates of volume change.

Data Analysis

We used the Post et al. (1971) tabular data to assess atypical changes in size of individual glaciers and as a check on our 1959 glacial extent. ArcGIS (10.2) was then used to determine glacial extent and area-altitude distributions (hypsometry) in 1959 and 2009. Glacial hypsometry at these times were obtained by overlaying glacier outlines on the 10 m DEM. ArcGIS was then used to determine glacierized area for the entire watershed in 50 m bins.

Glacial area in the bins was then multiplied by the surface melt from the monitored glaciers to estimate the magnitude and range in glacial runoff in 1959 and 2009 (Riedel et al. 2008). Summer melt on glaciers is strongly correlated with streamflow in this region (Stahl and Moore 2006, O'Neel et al. 2014) and decreases at higher elevations. We used the summer balance at 18 melt stakes on four glaciers distributed across the basin to define a regional summer melt vs. elevation curve. Surface mass balance measurements do not account for other sources of runoff from the glaciers, including summer precipitation and internal and basal ablation (Cogley et al. 2011). Melt season precipitation is typically minimal in this region, and surface melt is thought to account for roughly 90% of total glacier runoff (Mayo 1992).

An exponential curve provided the best fit to the melt data and was used to determine the amount of surface melt above and below the stakes. Glacier surface melt at the midpoint elevation of the 50 m bins was multiplied by the area of the corresponding bin to determine the volume of meltwater produced. The value for each bin was then summed to estimate total glacial runoff for the Skagit and three sub-basins.

Meltwater curves for extreme melt seasons (i.e. cool-wet vs. warm-dry) were combined with the 1959 and 2009 glacial hypsometry to assess the range of glacier surface meltwater production at each time. We chose melt data for 2003 and 2010 to represent extremes in glacial runoff because surface melt varied significantly between these years (∼ 1.5 m w.e. at a given point). The same melt curves were also used to determine melt from the 1959 hypsometry in the summers of 1963 and 1964 because only South Cascade Glacier was monitored in the basin before 1993. We assume the basin-wide glacial hypsometry did not change significantly from 1959 to 1964, and selected these years because they had extremes in summer balance at South Cascade Glacier (Josberger et al. 2007) and had similar summer positive degree-day values as 2003 and 2010.

The relative glacial surface meltwater contribution to total basin runoff during the summer melt season was based on a comparison of total streamflow measured at three USGS gaging stations on the Skagit River near Concrete and at the mouths of the Sauk and Baker rivers (Figure 2). The drainage area above the Concrete gage is about 7088 km², or 89% of the watershed. There are no glaciers and few tributaries that enter downstream of this gage.

Our ability to determine the relative contribution of snow, firn, and ice phases to total glacial surface meltwater production by this approach was limited by a lack of stake melt data for the mid-summer period. An estimate of volume change was used to make a rough estimate of the ‘fossil’ ice fraction contributing to summer glacial runoff (i.e. water lost from long term storage). Volume of glacial ice lost from 1959 to 2009 was estimated by multiplying the average cumulative annual mass balance by the glacier area in 10 year increments starting in 1979. The volume from each 10 year period was then summed to determine total volume loss. We assume there was minimal change in glacier volume from 1959–1979 because the cumulative surface annual mass balance of the five monitored glaciers was near zero (Riedel and Larrabee 2011). We also assume that the loss of glacier area was linear between decades since there were no other measures of area basin-wide between 1959 and 2009.

Errors

Historic and modem Skagit glacial surface meltwater volume estimates have several primary sources of measurement error: digitizing of glacier outlines, accuracy of the 10 m DEM, the hypsometry calculations with GIS, and the measure of surface melt at, above, and below the stakes. The uncertainty in glacier outlines includes random and systematic errors due to different human interpretation of the same imagery and the accuracy of the image. Our approach to calculating the area-altitude distribution of glaciers in 2009 using the 1959 DEM produces a slight over-estimate of glacier area at lower elevations because the DEM was not adjusted (lowered) for elevation changes associated with glacier recession and thinning. To assess this error we compared our 2009 hypsometry with data from the 2000 Shuttle Radar Topography Mission (SRTM) DEM.

Sources of error in the surface mass balance measurements were discussed in detail in Bidlake et al. (2010) and Riedel et al. (2008). Error in the measurement of summer melt (ablation) at the stakes was due to stake sinking and to estimates of the density and depth of snow and ice. Melt below and above the stake array adds additional uncertainty to our surface melt runoff estimates that was generally larger during warm summers. Assuming that the slope and position of the melt curves in the 1960s were similar to those from the 2000s introduces an unknown error that we qualitatively assessed by comparison of positive degree days (base 0°C) for the four years.

The estimate of total glacial ice volume lost from 1959–2009 has several sources of uncertainty, including errors in the measurement of glacier area change and cumulative annual mass balance. Surface mass balance measurements can significantly underestimate volume changes on temperate alpine glaciers if base maps used to integrate point (stake) measurements were not adjusted for surface elevation and glacier area changes. On South Cascade Glacier, this amounted to an estimate of cumulative annual mass balance that was 16% too positive in 27 years (Elsberg et al. 2001). Huss et al. (2009) note that for some glaciers surface mass balance measurements agree well with geodetic measurements. Errors in the estimates of total glacier area include those from the 1:24,000 scale topographic maps and the 2009 NAIP images. Maps used to reconstruct the mass balance of four glaciers from 1959 to 1992, and errors in the maps used to correct for bias in the surface measurements also add uncertainty.

Glacierized Area

Combined area of all Skagit Watershed glaciers larger than 0.1 km² in 1959 was 166.9 ± 8.3 km² (Post et al. 1971; Table 1). The accuracy of this inventory was reported as ‘excellent’ (5% error) on 379 of the 396 glaciers identified in the Skagit basin, while the 17 glaciers with ‘fair’ accuracy accounted for only 3.2 km² of the total glacierized area (< 2%). Our estimate based on the 1:24,000 scale map data from 1959 was 170.2 ± 8.5 km², but the difference with the Post et al. (1971) data was within the reported uncertainty for both data sets.

We estimate modem (2009) glacial cover in the Skagit basin at 138.21 ± 6.91 km², slightly less than identified by Dick (2013), but within the reported error of both surveys. The difference in digitized area was due to variable human interpretation of the same images that amount to an uncertainty of 2%. In a similar study in the Olympic Mountains, Riedel et al. (2015) estimated systematic error due to resolution of the NAIP images at 3%. Combined, systematic errors due to map resolution and debris/snow cover and errors due to human interpretation lead to an uncertainty in our 2009 combined glacier area of 5% (Table 1).

TABLE 1.

Change in combined glaciated area 1959–2009 for Skagit watershed and select sub-basins.

The number of glaciers in the Skagit basin larger than 0.1 km² declined from 396 in 1959 to 377 in 2009. Today, the Skagit watershed contains hanging and cirque glaciers, two small ice caps on Glacier Peak and Mount Baker, and an ice field southwest of North Klawatti Glacier (Figure 2). Boston Glacier, the largest in the watershed, shrunk from 7.00 to 6.11 km² since 1959. Nine glaciers that appeared to have grown based on comparison of the two data sets either actually got smaller (the three largest glaciers) or had positive area changes that were less than the 5% measurement error. In the case of the three largest glaciers, outlines on the 1959 DEM were inaccurate. For example, Quien Sabe Glacier (Figure 2) had an area of 0.67 km² on the 1959 map and 0.77 km² in 2009. When the glacier outlines were examined, it was clear that a large part of the glacier was missed on the 1959 map because it was obscured by snow. This interpretation was confirmed when noting that Post et al. (1971) listed the area of this glacier at 0.9 km² (i.e. it shrunk by 0.23 km² from 1959–2009).

Combined glacial extent in the Skagit basin decreased by 32.02 ± 14.42 km² (-19%) between 1959 and 2009. The uncertainty in this estimate was large because it includes the error from the 1959 and the 2009 data. There was a tendency toward less glacier area loss in more eastern parts of the Skagit basin (Table 1; Figure 2). Cascade River and the upper Skagit basin lost 17 ± 9%, whereas farther west the Sauk watershed lost 23 ± 9%. This pattern may reflect the fact that glaciers in the east were already restricted to high elevations above mean winter freezing lines in a colder continental climate. These glaciers also tend to have northern aspects that provide topographic shading.

Bolch et al. (2010) reported a decline of combined glacier area of 10% in the southern Coast Mountains of British Columbia, and 20% on Vancouver Island between 1985 and 2005. Riedel et al. (2015) observed a 34% decline in glacier area between 1980 and 2009 on the Olympic Peninsula. Although these other regional studies analyzed records that were 20–25 years shorter than in the Skagit, cumulative glacier annual mass balance was generally neutral to positive between 1959 and 1979 (Rasmussen and Conway 2001, Josberger et al. 2007). Therefore, glacier margins may not have changed significantly in the first part of our record, making it more comparable to other regional compilations.

O'Neal et al. (2015) used the same 2006 and 2009 NAIP images as we did to assess changes in the extent of 742 glaciers in the North Cascades, including those in the Skagit watershed. They found that 38% of the glaciers had gotten smaller, but that no statistically significant trend could be identified in the size of half the glaciers. Reliance on the 2006 NAIP images limited their ability to accurately digitize glacier margins in many cases because of snow that persisted on the landscape late into the summer. Our use of oblique photos and focus on the 2009 NAIP images allowed us to provide greater accuracy on many of the glaciers where they could not identify a change.

Figure 3.

Glacial hypsometry changes 1959–2009 in Skagit watershed. Dashed line denotes Shuttle Radar Topography Mission (SRTM) data obtained in 2000.

Glacial Hypsometry

Glacial area-altitude distributions in 1959 and 2009 are shown in Figure 3, with the Shuttle Radar Topography Mission (SRTM) data for comparison. The glacial hypsometry derived from the 2009 glacier outlines was very similar to the SRTM data, but our GIS approach tended to overestimate glacier area in the bins above 2100 m and underestimate it at lower elevation bins. We did not use the SRTM data because it was acquired in February, when several meters of snow covered North Cascade glaciers and because it was taken nine years earlier than the 2009 NAIP images, and missed several significantly negative mass balance years.

The largest loss of glacial cover occurred between 1750 and 2150 m elevation, where cirques comprise 15% of the high-elevation landscape in the North Cascades (Riedel and Probala 2005). Below 1650 m, the area loss by bin decreased steadily to a few tenths of a square kilometer. Factors that limit melt at the lowest elevations include valley orientation, debris-covered ice (mainly on Glacier Peak and Mount Baker), and shading of glacier termini in deep canyons. As expected, relatively little loss of glacier area occurred above 2200 m (Figure 3), which is slightly above the average equilibrium line altitude of glaciers monitored in the region (Osborn et al. 2012).

Glacial Mass Balance

The summer, winter, and annual mass balances of four glaciers that we monitor in the North Cascades are shown in Figure 4 (Riedel and Larrabee 2011). We do not include mass balance data collected after 2009 in our estimate of glacial runoff because of significant changes after the NAIP images were taken. Winter, summer, and annual surface mass balance measurement errors for the four NPS-monitored glaciers, averaged over the period of record, were ± 0.22, ± 0.29, and ± 0.35 m w.e., respectively (Riedel and Larrabee 2011). Box plots in Figure 5 show that differences in mass balance between the glaciers were close to these errors, but the data illustrate variability due to hypsometry, aspect, and distance from the Pacific Ocean. For example, Silver and Sandalee glaciers have less-negative average annual mass balances because they are located in a more continental climate, have north aspects, and are at higher elevations. As a result they get more of their precipitation as snow (higher winter accumulation) and have less summer melting. Noisy Glacier has a relatively high average annual mass balance for a glacier with a midpoint elevation of 1707 m because it gets additional winter accumulation from wind and avalanche deposition of snow from adjacent mountain slopes.

Figure 4.

Winter, summer and annual mass balance of four glaciers in and near Skagit watershed that are monitored by the NFS. See Figure 5 for average and range of values for individual glaciers, and Bidlake et al. (2010) for South Cascade Glacier data.

Figure 5.

Box plots of surface mass balance measurements on four glaciers in and near the Skagit basin from 1993 to 2013 (A-Noisy, B- Silver, C–N. Klawatti, D-Sandalee). Dark horizontal lines represent average balance, boxes represent the interquartile range (25 to 75th percentiles), and dashed lines are 5-95 percentile range of values.

All five monitored glaciers had significantly negative cumulative annual mass balance from 1959–2009 (Figure 6; Bidlake et al. 2010). From 1959 to 1979 the balance of all five glaciers was essentially neutral, and negative from 1979 to 2009. South Cascade Glacier had a cumulative loss of about -30 m w.e. due to its low elevation (Josberger et al. 2007), while North Klawatti lost less (-10 m w.e.) because of its higher elevation, large snow avalanche accumulation, and shaded aspect. North Klawatti Glacier had the least negative cumulative annual mass balance even though it had the most negative average annual mass balance because surface measurements were 10% more negative than the geodetically checked volume change from 1993–2006 (Figures 5 and 6; Riedel and Larrabee 2011). Potential reasons for the difference were that surface measurements miss large avalanche accumulation of snow on the south side of the glacier, and the maps used to correct the surface measurements for this glacier had large errors. Surface measurements on Silver and Sandalee glaciers were just the opposite, and were 16% and 14% too positive, respectively, from 1993–2006 (Riedel and Larrabee 2011). These values are comparable to the 17% overestimate of surface mass balance measured at South Cascade Glacier from 1970–1997 (Elsberg et al. 2001). Noisy Glacier required no adjustment to the mass balance measured on the surface. Net losses of volume from these glaciers is attributed to lower winter accumulation since 1976 (McCabe and Fountain 1995), lengthening of the melt season (Moore et al. 2009), and higher summer and winter air temperatures (Rasmussen and Conway 2001, Josberger et al. 2007).

Figure 6.

Cumulative annual surface mass balance of North Cascade glaciers 1959–2011. Mass balance reconstructed for four glaciers before 1998 by Granshaw (2001). See text for data sources.

Glacial Runoff

We refer to ‘glacial runoff’ specifically as the volume of meltwater derived from surface melt on all of the glaciers as measured at 18 sites on four glaciers (O'Neel et al. 2014). Our runoff estimates include snow melt on the glacier surface, as well as snow from previous winters (firn) and glacial ice. Most studies of glacial runoff focus on the ice melt component and assume that the contribution of snow and firn is minimal. While we do not focus on the relative contribution of each, we suggest that the snow and firn contribution should not be ignored since meltwater that arrives in June and early July is important during extended periods of low winter precipitation and/or high spring air temperatures. Further, even late in the melt season glacial runoff is typically a mix of snow, firn, and ice. Snow melt from a glacier's surface is also delayed in comparison to that of adjacent mountains slopes because glaciers benefit from cold-air pooling in cirques and create their own microclimates. Glaciers also absorb early season surface meltwater, and internal storage can delay its release by a month or more (Fountain and Tangborn 1985). In the case of firn, melt is delayed by a year or more.

Water years 1963 and 1964 were chosen to determine changes in glacial runoff from the 1959 glacial hypsometry, while 2003 and 2010 were selected to represent conditions a half-century later. These years were chosen because they represented extreme summer mass balance conditions and therefore glacial runoff (Figure 4; Josberger et al. 2007, Riedel and Larrabee 2011). Elevation vs. melt curves for 2003 and 2010 shown in Figure 7 were based on data from 18 melt stakes on four glaciers distributed across the basin, and indicate that summer 2003 had approximately one meter more melt at all elevations. Change in elevation was strongly correlated with the amount of surface melt in both low (r² = 0.85) and high melt years (r² = 0.87; Figure 7).

Figure 7.

Total summer (May–September) glacial surface melt at 18 stakes on four glaciers by elevation. Curves for cool-wet (2010—top) and warm-dry (2003) water years, and melt units are meters water equivalent.

The melt curves for 2003 and 2010 were assumed to be representative of the surface melt vs. elevation relationship for 1963 and 1964. This assumption was tested by comparing the sum of positive degree days (PDD base 0 °C) between the high glacial melt years (1963 and 2003) and the low melt years (1964 and 2010; Table 2). In both cases PDD values for high melt years were greater than low melt years, but summer 1964 had a PDD value 319 °C less than 2010. Summer 2010 was chosen to represent a low runoff year, yet the PDD value was only 86 °C less than the high runoff year 1963.

Fifty years ago glaciers produced substantial amounts of meltwater for the Skagit River, even though they covered only 2.4% of the watershed above Concrete (Table 2). In 1964, a water year with relatively cool, wet weather, the Skagit's glaciers produced approximately 0.440 ± 0.055 km³ of water during the summer melt season (Table 2). This represented about 6.7% of the total river discharge over the same period at the Concrete at 18 stakes on four top) and warm-dry equivalent. gage (Figure 2). In contrast, the previous summer of 1963 had warmer and drier weather that produced 0.742 ± 0.093 km³ of glacial surface meltwater. The 1964 glacial runoff estimate was high because we used a melt curve from a summer with 300 fewer positive degree days.

Errors from measurement of surface melt and hypsometry influenced our estimate of glacial runoff. Error in summer melt measurement averages ± 0.29 m w.e. on the four NPS glaciers from 1993–2011, or 6% of the melt at 2000 m elevation (Figure 7; Riedel and Larrabee 2011). Measurement of glacial area-altitude distribution in the basin includes error in the estimate of area (5%) and the DEM (1.5%; Fountain and Granshaw 2006), resulting in a total error in our runoff model of 12.5%. Melt below the elevation of the melt stake array adds additional uncertainty to our glacial surface melt runoff estimates (Figures 3 and 7). There was only 14.72 km² of glacierized area below 1700 m in 1959 (9% of total), and 8.84 km² (6%) in 2009. Significant glacial cover occurs above the highest melt stake near 2500 m, but relatively little melting occurs at these elevations (Figures 3 and 7). Post et al. (1971) estimated glaciers provided 0.8 km³ of melt water to the streams of the North Cascades region each summer, an area more than twice as large and with roughly twice the number of glaciers as the Skagit watershed. Drost and Lombard (1978), suggested that glaciers provided about 0.061 km³ to Skagit streamflow each year, based on a contribution of roughly 4.5 million m³ for every square kilometer of glaciated terrain. Our estimates (Table 2) for 1963 and 1964 are in the range of Post et al. (1971), but are about seven times higher than determined by Drost and Lombard (1978), who provided no explanation of how they determined glacier unit area runoff.

The loss of glaciers in the Skagit watershed in the last 50 years has cut deeply into summer base flow on the river. In the relatively warm, dry weather of 2003 glaciers contributed 0.559 ± 0.070 km³ of meltwater to total summer basin runoff at Concrete, a decrease of 24 ± 9% compared to 1963 (Table 2). These changes are significant at more than two times the error in our empirical model. Basin-wide recession of glaciers has also affected runoff in cooler, wetter years such as 2010, when glacier surface melt produced 0.333 ± 0.042 km³ of water, or 0.107 km³ less than 1964. Despite these losses, glaciers still contribute 6–12% of total summer streamflow to the Skagit River at Concrete (Figure 2). Regional studies have shown that glacial ice fraction of total basin summer runoff doubles in August and September when compared to the entire summer period (Tangborn and Fountain 1985, Bach 2003, Riedel et al. 2015). This means that in the Skagit the glacial surface meltwater fraction of total basin runoff was roughly 12–24% of the river flow in the last two months of the melt season.

TABLE 2.

Range and magnitude in the volume (km³) of the glacial surface melt contribution to total summer (May through September) streamflow for the Skagit River and its two largest tributaries. Four years shown represent strongly positive (cool and wet; 1964 and 2010) and negative (warm and dry; 1963 and 2003) glacial mass balance and weather years. Percent change based on comparison of high runoff years. See Figure 2 for locations of USGS streamflow gaging stations near Concrete (#12194000) and near the mouths of the Baker (#12193400) and Sauk (#12189500) rivers. Total streamflow in the upper Skagit was determined by subtracting the Baker and Sauk flows from that measured at the Concrete gaging station.

The change in summer glacial runoff across the watershed reflects the same pattern as the decline in glacial extent, although the difference was within the range of uncertainty. The pattern observed includes greater losses in the western Sauk basin (- 23 ± 9%) than in the more eastern upper Skagit (-16 ± 9%; Table 2). These changes are of the same magnitude as previous estimates in this region. In a study of several watersheds within and adjacent to the Skagit, Pelto (2008) identified a 25% decrease in glacial contribution to summer streamflow. Stahl and Moore (2006) found that negative trends in summer runoff from glacially fed rivers were widespread in southern British Columbia (Figure 2).

The loss of glacial surface meltwater to the Skagit River means that the overall variability in basin-wide runoff is increasing in early and late summer (Fountain and Tangborn, 1985, Stahl and Moore, 2006). Loss of glacial meltwater is important primarily later in the summer during periods of low river discharge caused by a lack of precipitation and snow-melt. Glacial meltwater from earlier in the melt season is likely offsetting a decline in snowpack (Mote 2003). These changes are expected to continue in the future, magnifying the importance of glacial sources of water (Elsner et al. 2010).

Intra-Basin Glacial Runoff

The influence of glaciers on summer streamflow varies considerably across the Skagit watershed depending on scale, glacial hypsometry, and distance inland from the Pacific Ocean source of moisture. More than half of the modem glacial surface melt was from the upper Skagit valley, including the Cascade River (Figure 2, Table 2). Glaciers in the Sauk valley provide another 20% and those in the Baker 24%. Suiattle River, White Chuck River, Baker River, Cascade River, and Thunder Creek remain as the primary glaciallyfed streams in the Skagit (Figure 8).

Thunder Creek in the upper Skagit has an extensive glacial cover of 13% and includes Fisher Creek and Thunder Creek HUC6 basins (hydrologic unit code; Figure 8; Seaber et al. 1987). This watershed has a favorable geographic position 100 km east of Puget Sound between the high-elevation Pacific and Skagit hydrologic crests (Figure 2). In this setting the basin receives the benefit of higher precipitation in a maritime-influenced highland climate, and the cooler accumulation season temperatures associated with the inland continental highland climate. We estimate Thunder Creek glacial surface melt runoff produced 0.107 km³ in summer 2003, or about 20% of the glacial total in Skagit River at Concrete. Despite the loss of 4.67 km² of glaciers in the last fifty years, Thunder Creek streamflow has changed less than six other large basins in the Puget Sound watershed (Pelto 2008) because the recession of glaciers has mitigated the loss of snow. Indeed, Frans (2015) found that summer runoff had increased on Thunder Creek from 1970 through 2009.

At a finer resolution, the Skagit has 75 sub-basins at the HUC6 scale (Figure 8). There are 31 HUC6 basins with no glaciers and 44 with glaciers, but 32 of those with glaciers have less than 3% glacial cover. Since 1959, Big Creek basin has lost all of its glaciers, although it is also one of the smallest at 55 km². Four other basins now have less than 0.1 km² of ice remaining and are likely to lose it in the next decade; all are located at the lower end of valleys (White Chuck River and Cascade River), or far to the east and on the leeward side of major hydrologic divides (upper Lightning Creek and Skagit headwater HUCs). The three HUC6 basins that lost the most glacial cover, and therefore had the largest reductions in glacial runoff, are Pass Creek (44%), Goodell Creek (41%) and the upper Whitechuck (33%). In terms of area lost, the eight basins that were the most glaciated in 1959 lost 15.33 km² by 2009, half of the total Skagit loss.

Volume Loss

The average cumulative annual mass balance for five Skagit watershed glaciers from 1959–2009 was - 20.35 ± 3.63 m w.e. and represents a layer of water lost from the surface of all the Skagit glaciers (Figure 6). Average cumulative annual mass balance was made significantly lower by the inclusion of data from South Cascade Glacier—the lowest elevation glacier monitored. We included it for our glacial surface meltwater estimate given the amount of ice lost in the Skagit basin at the altitude of 1650–2000 m (Figure 3). When the cumulative annual mass balance was multiplied by the decreasing area in 10-year intervals that were then summed, combined glacial volume loss between 1979 and 2009 was 3.01 ± 0.69 km³ w.e. The large error was due to uncertainty in the 1959 and 2009 maps (5%), the 1959–1992 reconstructed cumulative annual mass balance (10%; Granshaw, 2001), and in the maps used to correct the 1993–2009 surface mass balance estimates. Annual errors in the measurement of mass balance are compensating, and were determined by the square root of the sum of errors squared (Taylor 1997).

Granshaw and Fountain (2006) used scaling techniques based on glacier size/volume relationships to estimate glacial volume loss in North Cascades National Park at -0.8 ± 0.3 km³ from 1958–1998. The park is roughly one-third the size of the Skagit basin, contains 65 fewer glaciers, and encompasses the higher, eastern part of the basin where there are few large, low elevation glaciers such as those on Glacier Peak and Mount Baker (Figure 2). Considering these factors, the large error in both of our approaches, and the very negative annual mass balances from 2003–2009 (Figure 4), these estimates are comparable. Our volume loss estimates are within the range of those determined by Pelto (2008) in six watersheds within and near the Skagit basin. Volume losses in that study were determined by use of an average cumulative annual balance of -14 m w.e. between 1984 and 2004 and ranged from 20–40%.

Figure 8.

Fraction of glacial cover in Hydrologic Unit Code 6 scale sub-basins of the Skagit watershed in 2009.

Averaged over 30 years, the loss of glacial volume was 0.100 km³ each year. We assume that most of the volume lost was glacial ice (not snow), meaning that during the high melt year of 2003 roughly 18% of the total glacial surface meltwater volume reaching the Skagit at Concrete was from loss of ‘fossil water’ in the form of glacial ice. The ice melt fraction was higher in 2003, when bare ice was exposed earlier in the melt season than in 2010. At North Klawatti Glacier we have observed that the fraction of total summer surface melt due to glacial ice and not snow or firn varies from 6–33%. The ice melt component of glacial surface runoff is important because it reaches the Skagit River after mid-July when rainfall is minimal and snowpack is typically exhausted (Chennault 2001, Bach 2003).

Our estimate of a net loss of Skagit glacier volume from surface melt of 3.01 ± 0.69 km³ w.e. between 1959 and 2009 represents ∼ 100 years of fresh water supply for Skagit County at its current use rate of 0.029 km³ annually (7.7 billion gallons; Lane 2004). While undergoing significant changes during this period, the Skagit watershed continues to have more glacial cover than any other large basin in the US outside of Alaska (Fountain et al. 2007).