Study area

The southwestern Cordillera of Canada are comprised of numerous mountain chains including the southern Coast Mountains, the Columbia Mountains and the Canadian Rocky Mountains. Annual precipitation is delivered mainly in winter by low-pressure systems approaching from the southwest, producing strong west–east precipitation gradients and local distributions determined by windward and leeward slopes. Total annual precipitation averages 1230 mm at Whistler, British Columbia and 480 mm at Golden, British Columbia (Environment Canada, http://www.climate. weatheroffice.ec.gc.ca; Fig. 1) while mean net winter balances from 1965 to 1995 range between 1.75 m at Place Glacier and 1.19 m at Peyto Glacier (Reference DyurgerovDyurgerov, 2002). January (July) mean temperatures are −3.0°C (13.9°C) at Whistler (658 m a.s.l.) and −9.7°C (17.1°C) at Golden (780m a.s.l.), indicating that eastern sites will have a higher degree of continentality. Net glacier mass balances throughout the region have been predominantly negative since massbalance measurements began in 1965/66.

Fig. 1. Study area.

Mass-balance data

Mass-balance records at nine sites in southwest Canada (Fig. 1; Table 1) were obtained from published sources (Reference DyurgerovDyurgerov, 2002) and government reports (Reference Mokievsky-Zubok, Ommanney and PowerMokievsky-Zubok and others, 1985). Reported specific mass-balance data (b_{j ,z} ) for year j represent the winter (b _w) or summer (b _s) balance (in mm w.e.) at elevation z. For calculating positive degree-days, we assumed that reported specific mass-balance values represent the value for the midpoint of a given elevation band z. For convenience, summer balance is treated as a positive value.

Table 1. Mass-balance locations, with latitude (Lat), longitude (Long) and elevation range (z)

The period spanned by the mass-balance data includes substantial climatic variability. The period samples both cold and warm phase regimes of the Pacific Decadal Oscillation, as well as warm and cold phases of the El Niño–Southern Oscillation (Reference Moore and DemuthMoore and Demuth, 2001).

Initial analysis of mass-balance data indicated that some data filtering was necessary. Reported values that remained constant with elevation or followed fixed gradients (e.g. 500 mm w.e. (100 m)⁻¹) were removed from the analysis. These values were assumed to be estimated due to (1) loss of data associated with snow burial or meltout of ablation poles or (2) weather or logistical constraints preventing observation. We also removed observations where summer balance equaled zero (n = 69), as well as specific balance observations for Sykora Glacier that were duplicated from observations for Bridge Glacier (n = 55). The two glaciers were reported as one unit from 1977 to 1982 (Reference Mokievsky-Zubok, Ommanney and PowerMokievsky-Zubok and others, 1985). As a consequence of this data filtering, only individual specific mass-balance values were analyzed and no attempt was made to estimate mean specific summer balances.

Air-temperature data

Near-surface air temperatures were reconstructed for each glacier site using an inverse-distance-squared interpolation of federal and provincial climate data sources (Reference Stahl, Moore, Floyer, Asplin and McKendryStahl and others, 2006). To estimate average daily temperatures at the elevation bands of each site, we assumed a constant lapse rate of 6.0°C km⁻¹, which is consistent with previous melt-modeling studies in temperate locations (Reference MooreMoore, 1993; Reference Jóhannesson, Laumann and KennettJóhannesson and others, 1995; Reference VincentVincent, 2002).

Mass-balance sites used in this study are substantially higher in elevation than the climate stations (cf. Reference Stahl, Moore, Floyer, Asplin and McKendryStahl and others, 2006, fig. 1). However, low-elevation climate data are more likely to reflect variations in insolation and thermal characteristics of the regional air mass than climate data collected from within the glacier boundary layer, and thus may provide a more suitable positive degree-day (PDD) estimate (Reference Lang, Kraijenhoff and MollLang, 1986; Reference VincentVincent, 2002). In order that calculated melt factors were compatible with regional or gridded temperature fields, we did not consider glacier cooling effects.

The dates of mass-balance collection varied between sites and years, but were not reported for many observations. For the observations for which the measurement dates are known, the mean dates are 20 May and 28 September. To calculate PDDs, we assumed a nominal ablation season extending from 15 May to 30 September.

Derivation of melt factors

Separate melt factors for snow and ice can be fitted using an approach described below. The model was applied using PDD sums computed between 15 May and 30 September with a lapse rate of 6.0°C km⁻¹.

Below the glacier equilibrium-line altitude (ELA), specific summer balance is the sum of snow (M _s) and ice (M _i) melt (in mm), if melt from summer snowfall events is assumed to be negligible and ice melt occurs once the winter snow cover has been removed:

((1))

Seasonal snow and ice melt can be calculated as the product of melt factors (k _s and k _i) and the accumulated PDD for each surface (PDD_s and PDD_i):

((2))

((3))

((4))

where δ and ∊ are error terms.

Assuming that mid-summer snowfalls are small compared to the total winter balance, seasonal snowmelt (M _s) will be approximately equivalent to the amount of snow initially available for melting (or b _w):

((5))

Combining Equation (5) with Equation (2) yields

((6))

Inserting Equation (6) into Equation (4) and combining with Equations (1) and (2) yields an expression for the summer balance below the ELA:

((7))

which has the form of a multiple linear regression through the origin.

Above the ELA,

((8))

Equation (8) has the form of a simple linear regression through the origin. At the ELA,

((9))

Equations (7–9) suggest a piecewise linear regression model of the form of Equation (7) for b _w < k _s PDD, and of the form of Equation (8) for b _w ≥ k _s PDD. Melt factors for snow and ice were optimized for each site through non-linear least-squares regression using the ‘nlinfit’ function in Matlab. Summer balance b _s was the dependent variable and b _w and PDD were the independent variables. For each glacier, models were fitted using values for all years and elevation bands to generate one value for each of k _s and k _i. In addition, melt factors were estimated for each year for Peyto and Place Glaciers, which have the longest periods of record.

Sensitivity analyses of fitted melt factors

There are two sources of uncertainty in the calculation of PDD sums: (1) the length of the ablation season over which to calculate the sums (which should be based on the dates of mass-balance observations) and (2) the lapse rate used to adjust interpolated air temperature for elevation. To account for the former, four different ablation season lengths (ABL; Table 2) were used in the calculation of accumulated PDDs. Sensitivity of fitted melt factors to the method of temperature extrapolation was examined by calculating PDD totals using (1) a constant lapse rate of 6.5°Ckm⁻¹ and (2) mean monthly lapse rates which varied between 3.0 and 8.0°C km⁻¹, reflecting surface observations (Reference Stahl, Moore, Floyer, Asplin and McKendryStahl and others, 2006).

Table 2. Ablation season definitions used in calculating cumulative PDD

Another source of uncertainty in the fitted melt factors stems from the fact that the model does not distinguish between ice melt and firn melt when mass balances are negative. Firn is a winter accumulation of snow which persists through more than one melt season, giving it a lower albedo than fresh snow and a higher albedo than ice. While most degree-day models do not explicitly account for the presence or absence of firn, it is possible that the presence of firn might result in a lower fitted melt factor for ice. The piecewise model used in this study assumes that once the snow accumulations of the current winter are removed, ice melt begins. If firn were present, this would result in a reduced melt rate and thus a lower melt factor.

To evaluate the sensitivity of the fitted melt factors to the presence of firn, the surface type was estimated for each elevation band from cumulative and annual net balances (b _n). Observations where cumulative b _n was positive and annual b _n was negative were assumed to indicate that firn melting had occurred, and these data were flagged and removed from the statistical analysis. In addition, data from the first measurement year were removed if subsequent years demonstrated continuous firn cover. Piecewise regressions were then repeated at each site using the dataset filtered for firn melt.

Testing regional predictions of seasonal ablation

Specific summer balances at all sites were predicted using PDD _{j ,z} calculated from regional temperature interpolations and two sets of melt factors. The first set was derived from Place Glacier and was applied to assess the accuracy of summer balance predictions based on melt factors from a single long-term mass-balance site. The melt factors for Place Glacier were partly chosen because of the length of the mass-balance record, but also for its proximity to the majority of the sites studied. A second set of melt factors (k _s = 3.0 mm°C⁻¹ d⁻¹, k _i = 8.0 mm°C⁻¹ d⁻¹) was chosen to represent parameters assumed in several dynamic ice-sheet models (Reference Ritz, Fabre and LetréguillyRitz and others, 1996; Reference Marshall and ClarkeMarshall and Clarke, 1999; Reference Casal, Kutzbach and ThompsonCasal and others, 2004). These values are within the range of melt factors calculated for glacier sites using melt observations and in situ temperature measurements (Reference Braithwaite and ZhangBraithwaite and Zhang, 2000; Reference HockHock, 2003).

Seasonal snow- and ice-melt totals were estimated following Equations (2) and (4) by first using k _s and PDD _{j ,z} to deplete the observed winter balance b _w. Remaining positive degree-days were directed toward ice melt using k _i, and snow- and ice-melt totals (in mm w.e.) were summed to obtain an estimate of b _s (Equation (1)).

Static summer mass-balance sensitivity

The static summer mass-balance sensitivity (S_T ) to changes in temperature (Reference Oerlemans and FortuinOerlemans and Fortuin, 1992; Reference Braithwaite and ZhangBraithwaite and Zhang, 1999; Reference Braithwaite, Zhang and RaperBraithwaite and others, 2003; Reference De Woul and HockDe Woul and Hock, 2005) was determined by prescribing a daily 1 K increase, recalculating PDD totals and then re-estimating seasonal melt totals. The approach used here assumes no change in winter precipitation, and does not include the effect of summer snowfalls or a lengthened ablation season. Temperature sensitivity for each site was calculated as the mean difference between original and perturbed melt estimates, averaged over all data points.

Fitted melt factors

Fitted melt factors from the reference model run (15 May–30 September and a lapse rate of 6.0°C km⁻¹) range from 2.31 to 3.61 for k _s and 3.61 to 5.57 for k _i (Table 3). Except for Zavisha Glacier (where R ² is 0.37), model fits are reasonably strong with R ² up to 0.90. Calculated values of k _s and k _i are statistically different, with the exception of one outlying ice-melt factor derived from the weakest regression (Zavisha Glacier). Peyto Glacier has the strongest fit as well as the lowest k _s and the highest k _i.

Table 3. Melt factors for snow (k _s) and ice (k _i) and static mass-balance sensitivities (S_T ) to a 1 K temperature increase, calculated from the model run (15 May–30 September, 6.0°C km⁻¹)

Snow- and ice-melt factors are relatively consistent among sites, with mean (standard deviation) k _s =3.04 (0.38) and k _i = 4.59 (0.59). The mean values are similar to the optimized melt factors derived from a mass-balance and streamflow-calibrated hydrological model at Bridge Glacier (Reference Stahl, Moore, Shea, Hutchinson and CannonStahl and others, 2008), where k _s = 3.09 and k _i = 4.70 on 21 June. They also fall within the seasonal range of streamflow calibrated melt factors for the Lillooet River basin, located west of Place Glacier and north of Helm Glacier (Reference MooreMoore, 1993).

Melt factors calculated for Peyto Glacier (Fig. 2; Table 4) vary modestly from year to year, with low standard deviations and high correlation coefficients. Ice-melt factors exhibit higher variability between years than snowmelt factors. At Place Glacier, melt factors generally demonstrate a similar magnitude of interannual variability. However, the 1978 season shows k _s to be higher than k _i, suggesting a limit to this approach. In particular, fitting melt factors to individual years greatly reduces the degrees of freedom, and correspondingly increases the error in the fitted melt factors.

Fig. 2. Variation in annual melt factors calculated for Peyto (upper) and Place (lower) Glaciers.

Table 4. Statistical properties of melt factors (in mm°C⁻¹ d⁻¹) calculated annually for Peyto and Place Glaciers

Degree-day and melt factor sensitivity

Changing the ablation season length produced a non-linear response in accumulated positive degree-days with elevation for all lapse rate scenarios. Figure 3 illustrates two cases of this sensitivity for Bridge and Place Glaciers. At the higher elevations, PDD are nearly equal for all ablation season lengths. This is because mean daily temperatures are substantially greater than 0°C only on the warmest days (typically in June, July and August). At lower elevations, there is a larger spread between PDD totals calculated for different ablation seasons, but averaged over all elevations for Place Glacier there is only a 2% difference between the model run (15 May–30 September) and ABL3 (1 May–30 September). Extending the ablation season by 2months (1 April–30 October) increased accumulated PDD by an average of 4%, but only PDD accumulated between the dates of winter and summer mass-balance collection are of interest. Therefore, computed PDD are relatively robust to uncertainty in the dates of mass-balance observations. Examples of the sensitivity of PDD sums to lapse rate model are shown in Figure 4. For all sites and elevations, PDD calculated using varying monthly lapse rates are similar to those calculated using a fixed lapse rate of 6.0°Ckm⁻¹. For the examples shown, a fixed lapse rate of 6.5°C km⁻¹ decreased PDD accumulated between 15 May and 30 September by 13% for Bridge Glacier (1978) and by 15% for Place Glacier (1987) relative to the totals calculated using 6.0°C km⁻¹.

Fig. 3. Sensitivity of accumulated PDD to assumed ablation season length, calculated using a constant lapse rate of 6.0°C km⁻¹ for (a) Bridge Glacier, 1978 and (b) Place Glacier, 1987.

Fig. 4. Sensitivity of accumulated PDD to lapse rate model, calculated from 15 May to 30 September for (a) Bridge Glacier, 1978 and (b) Place Glacier, 1987.

Sensitivity analyses for the 12 PDD models (four ablation seasons and three lapse rate scenarios) for Peyto and Place Glacier melt factors are presented in Table 5. Melt factors are nearly equivalent for the variable lapse rate and the 6.0°C km⁻¹ fixed lapse rate, while a greater lapse rate results in higher fitted melt factors. Varying the time period over which PDD are calculated produces relatively small changes in calculated k _s, whereas k _i exhibits a stronger sensitivity to assumed ablation season length. However, even for k _i, the difference between the model run and the most likely range of ablation season lengths (ABL3, ABL4) is no more than 10% of the base case value. Accounting for the presence or absence of firn resulted in small and non-systematic changes to the fitted melt factors. Of the original 1029 mass-balance observations, 387 were removed due to the possible presence of firn and the melt factors were recalculated. At Peyto Glacier k _i increased by 4%, but the mean k _i for all sites increased by only 0.1% following the removal of mass-balance observations containing firn. Aerial photographs and the strongly negative cumulative mass balance of these sites suggest that firn extents are generally quite low in this region. However, melt factors for ice will be underestimated at sites where firn is more extensive, and in such cases mass-balance observations from the firn area should be excluded from the regression analysis.

Table 5. Results of sensitivity analysis for melt factors calculated from Peyto and Place Glacier mass-balance data (VAR: user-specified monthly lapse rate (see Reference Stahl, Moore, Floyer, Asplin and McKendryStahl and others, 2006, for details); 6.0 (6.5): constant lapse rate of 6.0 (6.5)°C km⁻¹). The reference model adopted for this study uses the 6.0°C km⁻¹ lapse rate and calculates PDD sums from 15 May to 30 September

Regional melt estimation

Melt factors assumed in several modeling studies (k _s = 3.0, ki = 8.0) generally over-predicted ablation. Melt factors fitted for Place Glacier (k _s = 2.59, k _i = 4.51) provided more accurate predictions of summer balance at most sites, although there was a tendency to under-predict at several sites including Bridge and Woolsey Glaciers (Fig. 5). Both sets of melt factors overestimated b _s for low summer ablation totals, and summer melt was poorly estimated at Helm and Zavisha Glaciers. The former was more accurately estimated using the higher ice-melt factor.

Fig. 5. Predicted vs observed summer balances for all sites. Summer balances modeled using k _i = 8.00, k _s = 3.00 are shown by crosses, and summer balances estimated with Place Glacier melt factors (k _s = 2.71, k _i = 4.69) are shown by circles.

Summer mass-balance sensitivity to a 1 K temperature increase

Calculated S_T for glaciers in western Canada range between −0.43 and −0.56 m w.e. a⁻¹ K⁻¹ (Table 3), which is within the range of values reported by other studies (Reference Oerlemans and FortuinOerlemans and Fortuin, 1992; Reference De Woul and HockDe Woul and Hock, 2005; Reference Braithwaite and RaperBraithwaite and Raper, 2007). These values are relatively low given the temperate location and moderate accumulation rates (average accumulation at the ELA of Place Glacier is 1.90m w.e.), but they do not include changes in the precipitation type due to increasing temperature. Maritime sites (Helm and Place Glaciers; Fig. 1) tend to have higher sensitivities than continental sites (Peyto Glacier), which is consistent with previous studies.