A preliminary assessment of glacier melt-model parameter sensitivity and transferability in a dry subarctic environment

Efforts to project the long-term melt of mountain glaciers and ice-caps require that melt models developed and calibrated for well studied locations be transferable over large regions. Here we assess the sensitivity and transferability of parameters within several commonly used melt models for two proximal sites in a dry subarctic environment of northwestern Canada. The models range in complexity from a classical degree-day model to a simplified energybalance model. Parameter sensitivity is first evaluated by tuning the melt models to the output of an energy balance model forced with idealized inputs. This exercise allows us to explore parameter sensitivity both to glacier geometric attributes and surface characteristics, as well as to meteorological conditions. We then investigate the effect of model tuning with different statistics, including a weighted coefficient of determination ( wR2), the Nash-Sutcliffe efficiency criterion (E), mean absolute error (MAE) and root mean squared error (RMSE). Finally we examine model parameter transferability between two neighbouring glaciers over two melt seasons using mass balance data collected in the St. Elias Mountains of the southwest Yukon. The temperature-index model parameters appear generally sensitive to glacier aspect, mean surface elevation, albedo, wind speed, mean annual temperature and temperature lapse rate. The simplified energy balance model parameters are sensitive primarily to snow albedo. Model tuning withE, MAE and RMSE produces similar, or in some cases identical, parameter values. In twelve tests of spatial and/or temporal parameter transferability, the results with the lowest RMSE values with respect to ablation stake measurements were achieved twice with a Correspondence to: G. E. Flowers (gflowers@sfu.ca) classical temperature-index (degree-day) model, three times with a temperature-index model in which the melt parameter is a function of potential radiation, and seven times with a simplified energy-balance model. A full energy-balance model produced better results than the other models in nine of twelve cases, though the tuning of this model differs from that of the others.


Introduction
Climate warming is expected to reduce the extent of Earth's mountain glaciers and ice caps during the 21st century, raising eustatic sea level and diminishing fresh water resources (e.g.Lemke et al., 2007).In the past decade there have been attempts to project the magnitude of glacier loss using melt models applied over large regions or even globally (e.g. de Woul and Hock, 2005;Oerlemans et al., 2005;Raper and Braithwaite, 2006;Hirabayashi et al., 2010;Radić and Hock, 2011).Such studies have produced a wide range of projected contributions of mountain glaciers and ice caps to 21st century sea-level, from 4 cm Sea Level Equivalent (SLE) (Raper and Braithwaite, 2006) to 36 cm SLE (Bahr et al., 2009).Contributing to this range are uncertainties in the total volume of glaciers and ice caps (e.g.Raper and Braithwaite, 2005;Meier et al., 2007), variation in the output from different climate models (Randall et al., 2007;Radić and Hock, 2011) and (the focus of this manuscript) assumptions made in applying melt models outside of the domain where they were developed and tested.
Conservation of energy is the physical law that controls the melting of snow and ice (Oke, 1987).For a given volume this law can be written as: A. H. MacDougall et al.: Glacier melt model parameter sensitivity and transferability where Q in is incoming energy flux, Q out is outgoing energy flux, and Q store is a storage term that also takes into account energy production within the system.Equation ( 1) is here referred to as the "true energy balance".For the purposes of studying and modelling energy balance, the incoming and outgoing fluxes are broken into energy balance terms.
The energy balance terms source energy exchange to specific physical processes operating on and within the system.The process of breaking the energy balance into terms and the description of the physical processes underlying these terms is only ever approximate but is nonetheless a useful way of conceptualizing the melt of ice and snow (Oke, 1987).Glacier melt models have been broadly divided into temperature-index models, which correlate melt to air temperature, and energy-balance models, which use energybalance theory to solve for the energy available to melt snow or ice (e.g.Hock, 2003Hock, , 2005)).Both temperature-index and energy-balance models attempt to approximate the true energy balance, the former through empirical parameterization at the level of melt and the latter through process-based modelling of individual terms in the energy-balance.If an energybalance model gives a good approximation of the true energy balance, then this model can be used to study the nature, sensitivity and transferability of simpler empirical models.Temperature-index models have been preferred for large scale application due to their low data requirements and overall good performance (e.g.Huss et al., 2008).A review of melt factors (the controlling parameters in temperature-index models) from difference studies, however, reveals large variations between regions with no obvious climatically based pattern (Hock, 2003).This result hints at the issue of parameter transferability and suggests that caution should be used when applying these models globally.
Model transferability is a measure of the generality of a model.A model is considered transferable if it is able to produce realistic results outside of the domain for which it was developed and tested (Takle et al., 2007).Model transferability is closely related to the concept of model sensitivity, which is a measure of how variation in model output can be attributed to variation in model input (Saltelli et al., 2004).Transferability is a function of the sensitivity of the model to changes in environmental conditions, and the difference in environmental conditions between the domain in which the model is developed versus the domain in which it is implemented.Model sensitivity arises from multiple sources, including the structure of the model, the population and representation of processes in the model and the values of model parameters (Saltelli et al., 2004).The transferability of a given glacier melt-model can be broken into two categories: parameter transferability (the applicability of model parameters to a time or location other than those where they were derived or measured) and meteorological transferability (the applicability of meteorological conditions measured at one location to another).Meteorological transferability is a well studied field (e.g.Wilks, 2006) and we will therefore focus exclusively on model parameter transferability.Various metrics can be used to quantify model parameter transferability, the simplest being the parameter values themselves.Identical models that have identical parameter values will be fully transferable.
Two recent studies have examined the transferability of glacier melt-model parameters, one for the mountains of southwestern Canada (Shea et al., 2009) and the other for the Swiss-Italian Alps (Carenzo et al., 2009).Shea et al. (2009) examined the stability of melt factors for a temperature-index model applied to nine glaciers in the southwestern Canadian Cordillera.The melt factors were found to be highly consistent between the glaciers despite the varying maritime to continental climate settings.Carenzo et al. (2009) examined the transferability of several variations of the enhanced temperature-index model of Pellicciotti et al. (2005).The model parameters were found to be highly transferable between five point locations on the same glacier in one ablation season, between three ablation seasons at one location on the same glacier, and for two additional point locations on two other glaciers.
Here we expand on previous work addressing glacier meltmodel transferability in the Donjek Range of the St. Elias Mountains (MacDougall and Flowers, 2011) by extending the analysis of model parameter transferability to several commonly used glacier melt models.This analysis is done in two stages: (1) the sensitivity of model parameters to different environmental conditions is explored by tuning the models to the output of an energy-balance model forced with idealized inputs; (2) the transferability of melt-model parameters between two glaciers and over two melt seasons is examined by tuning the models to mass balance data collected in the subarctic study region.

Study site
In this study we target two individual valley glaciers in the Donjek Range of the St. Elias Mountains, southwest Yukon, Canada (Fig. 1).The St. Elias Mountains are extensively glacierized (Arendt et al., 2008) and have contributed significantly to sea level during the latter half of the 20th century (e.g.Kaser et al., 2006;Berthier et al., 2010).The Donjek Range represents a transitional region between the ice-free foothills to the northeast and the contiguous icefields to the southwest.The climate in this region can be characterized as subarctic, due to the strong orographic blocking of the St. Elias Mountains (Marcus and Ragle, 1970).
The glaciers of interest lie between 60 • 47 N and 60 • 57 N, and 139 • 05 W and 139 • 13 W.Although they are unnamed, we refer to them here as "South" and "North" Glaciers, indicating both the respective sides of the local range crest on which they are situated and their dominant aspects.South Glacier has an area of 5.3 km 2 and spans an elevation range of 1970-2960 m above sea level (a.s.l.).It is known to be polythermal based on numerical modelling (De Paoli and Flowers, 2009), several measured englacial temperature profiles and extensive radar survey (unpublished data, Simon Fraser University Glaciology Group).The glacier is known to be surge-type (Johnson and Kasper, 1992), and is thought to be undergoing a slow surge at present (De Paoli and Flowers, 2009).North Glacier has an area of 6.9 km 2 and ranges from 1890-3100 m a.s.l. in elevation.It also has a polythermal structure based on ice temperature measurements and radar data, but is not known to surge.The present equilibrium line altitude (ELA) for both glaciers is approximately 2550 m a.s.l.(Wheler, 2009).

Field data collection and processing
Automatic Weather Stations (AWS) were deployed at ∼2300 m a.s.l. in the ablation zones of North and South Glaciers from 2007 until 2009, with a full complement of instruments deployed in May 2008.The AWSs are instrumented to measure air temperature, barometric pressure, relative humidity, wind speed, rainfall, net all-wave radiation, and incoming and outgoing shortwave radiation (Table 1).Instruments other than the rain gauge are installed at a nominal height of 2 m on a tripod that sits on the ice surface.Instruments thus maintain a relatively constant height above the surface.Except for barometric pressure and rainfall (which are measured every half-hour) the measurements are taken at five minute intervals.These data are gap-filled using linear interpolation, however no gap is longer than 30 min for deployed and undamaged instruments.The gap-filled data are averaged to hourly values (except rainfall where hourly totals are used).An ultra-sonic depth gauge (USDG) is located several meters from each AWS.The USDGs measure the distance to the surface and are used to estimate snowfall during the summer season.A snowfall event is interpreted to have occurred if the daily-average distance to the surface decreases.New snow is assigned a density of 200 kg m −3 based on field measurements.Table 2 summarizes the mean meteorological conditions during the 2008 and 2009 melt seasons.We pragmatically refer to the summer season for each year as the time interval during which the albedometers are deployed at both sites.For 2008 this period is 6 May to 14 September, while for 2009 it is 8 May to 2 August.
An array of 17-18 ablation stakes is maintained on each glacier (Fig. 1).The height of the stakes above the surface and the surface density are measured at the beginning and end of the summer season on North Glacier and at weekly to monthly intervals on South Glacier during the melt season.These height and density measurements are converted into mass balance estimates using the methods of Østrem www.the-cryosphere.net/5/1011/2011/The Cryosphere, 5, 1011-1028, 2011  (1991).Ice is assumed to have a density of 900 kg m −3 .Snow pits are excavated to assess the density structure of the snowpack in the ablation and accumulation zones of each glacier in May of each year (before the onset of melt).The snow pit data and the initial snow depth at the ablation stake locations are used to estimate the winter accumulation.For South Glacier, this is done by regressing the water-equivalent snow depth against surface slope and elevation.For North Glacier, where surface slope variations are much less than on South Glacier, the snow depth relationship is adequately captured by regression on elevation alone.Summer snowfall events detected by the USDGs are extrapolated from the USDG locations to the rest of the glacier using a precipitation lapse rate (Table 3) and assuming precipitation falls as snow when air temperature is below a threshold of 1 • C (Jóhannesson et al., 1995).The precipitation lapse rate is based on field measurements (see Wheler, 2009).The mean firn line elevation of each glacier is estimated based on field observations during ablation stake surveys.Temperature lapse rates are measured within the glacier boundary layer by Onset HOBO ™ microloggers deployed on a subset of the ablation stakes (Fig. 1).The temperature values recorded by each logger are averaged for the whole summer season and linear regression used to compute a lapse rate.Our data permit this to be done for both glaciers in 2008, but only for South Glacier in 2009.

Melt models
A large number of formulations for empirically based glacier melt models have been proposed in the literature (see Hock, 2003 for review).Here we focus on some of the more commonly used formulations that can be employed in a spatially distributed fashion and that have previously been used in studies of melt model transferability or long-term projections of glacier mass balance (Shea et al., 2009;Carenzo et al., 2009;Hock et al., 2007).The four glacier melt models we include are: (1) the classical temperature-index model (CTIM) (Braun et al., 1993), (2) the temperature-index model of In addition to the models above, we outline the energy balance model (EBM) of MacDougall and Flowers (2011), used variously in this study for tuning the other models and for the purposes of model-output comparison as detailed below.

Classical temperature-index model (CTIM)
The CTIM correlates air temperature to melt with an empirical degree-day factor.As in most other implementations of this model, separate degree-day factors are used for ice and snow (e.g.Braun et al., 1993).The transferability of this model (via its degree-day factors) has previously been examined by Shea et al. (2009) for glaciers in southwestern Canada.The model takes the form: where M is melt rate and DDF snow/ice is the degree-day factor for snow or ice.The CTIM is driven with air temperature T a and precipitation.We adjust T a for elevation by prescribing a constant temperature lapse rate.The treatment of precipitation varies by experiment and is explained in the relevant sections below.Firn is treated as snow in the model by assigning an arbitrarily deep snow depth above the firn line.
The models described below treat air temperature, snowfall, and firn identically to the way they are treated in the CTIM.

Temperature-index model of Hock (1999) (HTIM)
The HTIM is an extension of the temperature-index method, where the degree-day factor is parameterized as a linear function of potential shortwave radiation.This model has been widely used and exhibits significant improvements in predictive capability over the classical temperature-index approach, The Cryosphere, 5, 1011-1028, 2011 with a minimal increase in data requirements (e.g.Hock, 1999;Huss et al., 2008).The model takes the form: where MF is a temperature melt factor, r snow/ice is the radiation melt factor for snow or ice, and I p is the potential direct shortwave radiation.I p varies in time and space due to the combined effects of the position of the sun, surface slope, surface aspect, and shading from surrounding topography.I p is calculated for each grid point using solar geometry and digital elevation models (DEMs) of the glacier and surrounding terrain, and by assuming a constant diffuse fraction of radiation.It has been noted that the form of this model, in multiplying temperature and shortwave radiation together, is physically problematic (Greuell and Genthon, 2004).Caution is therefore advised when interpreting the model results.

Temperature-index model of Pellicciotti et al. (2005) (PTIM)
The PTIM uses an arithmetic combination of terms representing the contributions of incoming shortwave radiation and air temperature to melt.The model was developed by Pellicciotti et al. (2005) and has been applied at the point scale for glaciers in the Swiss and Italian Alps and in the semi-arid central Andes of Chile (Carenzo et al., 2009;Pellicciotti et al., 2008).The model takes the form: where α is albedo, S in is incoming shortwave radiation, T F is the temperature melt factor, and S RF is a shortwave radiation melt factor.The treatment of the incoming shortwave radiation (S in ) and albedo (α) are described in the context of the relevant experiment.

Simplified energy-balance model (SEBM)
The simplified energy-balance model of Oerlemans (2001) takes the form: where Q M is the energy available to melt snow or ice, L f = 3.34 × 10 5 J kg −1 is the latent heat of fusion, ρ w = 1000 kg m −3 is the density of water and C 0 and C 1 are empirical factors that together take into account net longwave radiation and the turbulent heat fluxes.Note that the model does not employ a temperature threshold and therefore that the temperature term is negative when air temperature is negative.As with the PTIM the treatment of incoming shortwave radiation (S in ) and albedo (α) varies by experiment and is described below.

Energy balance model (EBM)
The energy balance at a glacier surface is often expressed as: where the heat flux due to rain is neglected, and S in is incoming shortwave radiation, α is the surface albedo and L in and L out are the incoming and outgoing longwave radiation, respectively.Q H is the heat flux due to the difference in temperature between the glacier and the atmosphere (sensible heat flux), while Q L is the heat flux to or from the glacier via material phase change at the surface (latent heat flux).Q g is the energy released or absorbed by the subsurface when the snow or ice changes temperature.
The Stefan-Boltzmann relationship and surface temperature (T s ) are used to compute outgoing longwave radiation (L out ).A simplified subsurface scheme is used to calculate the surface temperature and the subsurface heat flux (Q g ).This scheme forces the subsurface heat flux into a thin layer when the residual of the energy balance is negative: where t is the model time-step, ρ s is the surface density, the specific heat capacity of ice (c s ) is equal to 2110 J kg −1 K −1 , and d s = 0.1 m is the prescribed thickness of the subsurface layer.This scheme is a compromise between a more complicated multi-layer subsurface model and simpler iterative approximations or the assumption of a constant surface temperature (e.g.Wheler and Flowers, 2011).It also allows for temporary heat storage in the subsurface with minimal data requirements (Wheler and Flowers, 2011).The treatment of the remainder of the radiative balance differs in the real and idealized implementations of the model and is therefore described later.
The bulk aerodynamic approach is used to calculate sensible (Q H ) and latent (Q L ) heat fluxes (e.g.Anderson et al., 2010;Anslow et al., 2008;Brock et al., 2000;Hock and Holmgren, 2005).The fluxes are described as: where ρ a is the density of air, c p = 1004 J K −1 kg −1 is the heat capacity of air, k = 0.4 is the von Kármán constant, u z is wind speed, z is the measurement height above the surface, z o is the aerodynamic roughness length of the surface, z oT is the roughness length for temperature, L is the Obukhov length, M,H are stability constants, e z is the vapour pressure at height z, e s is vapour pressure at the surface, L v is the latent heat of vapourization, p c is the atmospheric pressure and z oe is the roughness length for humidity.An iterative loop is needed to solve for the turbulent fluxes because L is a function of Q H .The aerodynamic roughness-length (z o ) used in the bulk aerodynamic approach is either modelled separately or prescribed as a constant.The roughness lengths for humidity and temperature are taken to be two orders of magnitude smaller than that for momentum, following Hock and Holmgren (2005).

Implementation with idealized inputs (IEBM)
To investigate melt-model parameter sensitivity to glacier geometry, surface conditions and meteorological variables, we tune the melt models under consideration to the output of the EBM forced with idealized inputs.We refer to this implementation of the energy balance model as the "IEBM".
In the IEBM, incoming shortwave radiation (S in ) is computed from top-of-the-atmosphere radiation multiplied by a factor B x that accounts for atmospheric absorption and reflectance.Top-of-the-atmosphere radiation is computed using well known equations for solar geometry (Oke, 1987).Shortwave radiation is broken into direct and diffuse components using a constant partitioning ratio D f .Diffuse radiation is imposed on all grid cells equally, while direct radiation is adjusted according to the slope and aspect of each grid cell.Ice and snow surface albedos, α ice and α snow , are taken as constant in the IEBM and based on averages of these quantities from the Donjek Range study sites (Table 4).In the calculation of the turbulent heat fluxes Q H and Q L , roughness length z o is also taken as a constant (see Table 4).
Incoming longwave radiation (L in ) is computed using the parameterization of Greuell and Knap (1997) which takes the form: where sky is the emissivity of the sky, σ = 5.67 × 10 −8 J s −1 m −2 K −4 is the Stefan-Boltzmann constant, cs is the clear sky emissivity, n is the cloud fraction, p = 2 is an exponent, cl = 0.983 is the cloud emissivity and b = 0.433 is a constant.The values of these parameters are taken from Greuell and Knap (1997) without modification.This parameterization has previously been used by Klok and Oerlemans (2002).
The only time-varying meteorological inputs to the IEBM are incoming shortwave radiation (see above) and temperature.All other meteorological variables, including wind speed, barometric pressure, cloud fraction and relative humidity, are held at constant values for ease of interpretation (see Table 4).Temperature is represented with two superimposed cosine curves as follows: T a (t)=µ 1 cos 2π 365.24 (t+µ 2 ) +µ 3 cos(2π (t+µ 4 ))+µ 5 , (14 where t is time in days and µ 1:5 are, respectively, annual temperature amplitude and phase, daily temperature amplitude and phase, and annual mean temperature. The IEBM employs an idealized glacier geometry wherein the ice surface is described by a linear equation of position and elevation, with inputs of maximum (Z max ) and minimum (Z min ) elevations, glacier length (L g ) and aspect (A).The glacier surface is one dimensional.Glacier slope is constant and computed from the inputs above.Initial snow depth is prescribed as a function of surface slope (C s ) and elevation (C z ) with intercept C I .Precipitation and firn are neglected in the IEBM.When the snowpack is removed via ablation, ice is assumed to be exposed regardless of elevation.

Implementation with real data (REBM)
For the purposes of comparison, we highlight previously published results of the energy balance model forced with real data (MacDougall and Flowers, 2011).We designate this implementation of the model "REBM".The comparison between results of the REBM and the other melt models should be interpreted with caution, as the REBM is not tuned in the same fashion as the other models: rather than tuning In the REBM, shortwave radiation is broken into its direct and diffuse components using the method of Collares-Pereira and Rabl (1979) and Hock and Holmgren (2005).Direct shortwave radiation is only incident on the fraction of the glacier unshaded by surrounded terrain, while diffuse radiation is assumed to originate from all parts of the sky equally and is applied to all grid cells.If the AWS is shaded by surrounding topography, all measured incoming shortwave radiation is diffuse (Hock and Holmgren, 2005); in this situation the ratio of direct to diffuse shortwave radiation from the most recent time where the AWS was unshaded is used to approximate direct shortwave radiation at unshaded grid cells.In practice, this only occurs when the sun is close to the astronomical horizon when shortwave radiation is weakest.Skyview fraction and topographic shading are computed using DEMs of the terrain surrounding the glaciers and (for topographic shading) the traverse of the sun through the celestial hemisphere (e.g Oke, 1987).
The REBM uses the albedo parameterization of Hock and Holmgren (2005): where α t is the albedo at time t, n d is the time since the previous snowfall in days, P s is the snowfall rate and a 1:4 are constants that are found by tuning the model to a measured albedo record (Hock and Holmgren, 2005).Ice is assumed to have a constant albedo that is taken as the mean value of ice albedo measured at the AWS location (α ice ) (e.g.Hock and Holmgren, 2005;Oerlemans and Knap, 1998).
Incoming longwave radiation in the REBM is inferred from measurements of shortwave and net all-wave radiation, along with modelled surface temperature from a multilayer subsurface heat-flux model applied only at the AWS location (Wheler, 2009).Turbulent heat fluxes Q H and Q L , along with the subsurface flux Q g , are simulated in an identical fashion to the IEBM, except that the evolution of the aerodynamic roughness length of snow is parameterized as in Brock et al. (2006).The aerodynamic roughness length of ice is taken as constant and equal to the mean of the measured roughness lengths of ice.See MacDougall (2010) and MacDougall and Flowers (2011) for details, including model validation and sensitivity analysis.
In contrast to the IEBM, the REBM uses measured timeseries of incoming shortwave radiation, temperature, relative humidity, wind speed, barometric pressure and precipitation.It employs the real glacier DEMs and corresponding calculations of topographic shading and sky-view fraction, as well as field-based estimates of temperature and precipitation lapse rates.For North Glacier in 2009, where the data required to estimate temperature lapse rate are unavailable, the value for 2008 is used.

Optimization statistic
An empirical model must be tuned to a data set that the model is capable of simulating in order to obtain model parameter values (Krause et al., 2005).This is often accomplished by running the model for a spectrum of parameter values and comparing the model output to the data using some statistic.In this experiment we examine the effect of the choice of optimization statistic on tuned parameter values for each of the melt models under consideration.The four statistics examined are: (1) the weighted coefficient of determination (wR 2 ), where the weight w can take on values between 0 and 1: (2) the Nash-Sutcliffe efficiency criterion: (3) the mean absolute error: and (4) the root mean square error: where M s i is the simulated ablation and M r i is the reference or measured ablation at location i, n b is the number of locations at which ablation values are available and g is the slope of the regression upon which R 2 is based.R 2 also produces an intercept that ideally should be zero for a 1 : 1 correlation.Note that both E and RMSE contain the sum of the squared differences.
Each of the melt models is tuned to the cumulative ablation as simulated by the IEBM with the idealized inputs intended to represent Donjek Range glaciers and their environment (Table 4, column labelled "DRG").Our tuning method simply involves discretizing the parameter space of a model and running the model with all possible parameter combinations within a plausible range of values for each parameter.The best match between simulated and reference ablation is found by maximizing wR 2 or E (maximum value is 1), or minimizing MAE or RMSE (minimum value is 0).This procedure is repeated for glacier aspects aligned with the four cardinal directions.In all experiments where the melt models are tuned to IEBM output, the SEBM uses S in and static albedos identical to that of the IEBM.This gives the SEBM an advantage over the other models.

Parameter sensitivity
Glacier geometry, surface characteristics and meteorological conditions all affect the surface energy balance.Here we examine the effects on melt model parameters of glacier aspect, surface slope, mean, minimum and maximum elevations, snow and ice albedo, snow depth at the onset of the melt season (winter balance), wind speed, annual mean temperature and lapse rate.We do this by altering each parameter independently within the IEBM and then tuning the melt models to the cumulative melt predicted by the IEBM using RMSE as the optimization statistic.The range of values chosen for each characteristic encompasses the differences between the two study glaciers and two melt seasons under consideration (Table 2), and is extended to capture the plausible range of values for the region.Certain geometric attributes cannot be manipulated without changing others.Slope is manipulated jointly with glacier length to maintain a constant elevation range.Minimum and maximum elevation are changed independently which changes the mean elevation in these tests.Where the mean elevation is manipulated, maximum and minimum elevations also are changed, but slope and glacier length are held constant.For the mean elevation test, the intercept of the initial snow-depth relationship is changed in tandem with mean elevation, in order to prevent unrealistic snow depth values.The barometric pressure (taken as a constant in the IEBM) is changed in the mean elevation test in accordance to the United States standard atmosphere of 1976 (National Aeronautics and Space Administration et al., 1976).Two tests are conducted using lapse rates: one in which the melt models use a lapse rate identical to that in the IEBM, and one in which the melt models use a fixed lapse rate of −6.5 K km −1 .For the purposes of graphical comparison, the parameters from each of the melt models are converted to a common set of units (m w.e.).This is accomplished by multiplying each parameter by its index of melt energy (e.g.positive degree days for the CTIM) calculated on an unshaded horizontal reference surface at a point 2300 m a.s.l. in the IEBM.

Parameter transferability
Parameter transferability is investigated using real data, with the melt models being tuned to cumulative ablation measured at stake locations rather than the output of an EBM.For each stake only the ablation between the first and the last measurement of the season are used in order to avoid autocorrelation.The models are driven with air temperature measured at the AWS locations and are implemented in a distributed fashion using the glacier surface DEMs.Shortwave radiation in the SEBM is treated identically as in the REBM: S in is measured at the AWS locations and albedo model parameters for each glacier and year are found by minimizing RMSE between simulated and measured albedo values.A constant firn-line elevation of 2450 m a.s.l. is prescribed for all models, above which the snowpack is assumed to be arbitrarily deep.Summer snowfall is extrapolated as described in Sect.3.1.Melt-model tuning is initially performed separately for each glacier (North, South) and year (2008,2009), yielding four sets of parameters for each model that define the local "control" runs.One additional "master" control run is performed for each model using the full complement of data from both glaciers and both melt seasons.The master control run has a similar design to some of the earlier work on glacier melt model transferability (e.g.Braithwaite, 1995).The control runs collectively serve as references against which to evaluate the results of the transferability tests.
We assess melt-model parameter transferability in time, in space and in space and time together for each model (Fig. 2).In each test, the parameter values derived for one glacier and year are used in place of those locally derived for the other glacier and/or year.In transferring parameters between glaciers and/or years we aim to test the hypothesis that melt can be accurately modelled with parameters derived from other sites or derived locally in other years (as has been shown by Shea et al. (2009) in southwestern Canada and by Carenzo et al. (2009) in the Swiss-Italian Alps).We use the RMSE between the simulated and measured cumulative ablation at the stake locations to evaluate the success of the control runs and parameter transfer tests.

Optimization statistic
In this experiment, two of the four statistics (E and RMSE) produce identical parameter values for each model with our tuning method (Table 5).Both E and RMSE contain the www.the-cryosphere.net/5/1011/2011/The Cryosphere, 5, 1011-1028, 2011 w.e.µm h sum of the squared differences as a key part of their definitions.Using MAE produces results that are very similar to those using E and RMSE, varying by one interval of parameter discretization at most for the CTIM, PTIM and SEBM, and two intervals at most for the HTIM.However, using MAE leads to MF = 0 for the HTIM in every test in Table 5.The HTIM is known to exhibit equifinality (or nonuniqueness) in model parameters, such that multiple parameter sets can yield equally good model performance (Carenzo et al., 2009).This may explain the more variable response of the HTIM and suggests caution in interpreting the results of this model.Bearing this in mind, Table 5 reveals little difference in the values of optimized parameters when tuning with E, RMSE or MAE and therefore little difference in the ablation amounts these models would predict when tuned with these statistics.Employing the conditions used to convert the melt parameter units to m w.e.(see Sect.  4).
MAE results in less than a 5 % variation in simulated ablation.
The results for wR 2 differ greatly from those of the other statistics.This disparity is related to the fact that R 2 only quantifies the variability of observed and predicted data.Simulations can have high R 2 values but systematically underestimate or overestimate ablation (Krause et al., 2005).We have attempted to compensate for this flaw by weighting R 2 by the slope of the regression between observed and predicted values (Eq.17), but this has not eliminated the systematic bias as the intercept of the relationship can be nonzero.

PTIM tuning and the effect of climate setting
For synthetic glaciers oriented in three of the cardinal directions, when the PTIM is tuned to IEBM-derived cumu-lative ablation for Donjek Range climate conditions, the optimal value of the solar radiation factor (S RF ) is found to be zero when tuning with E, RMSE or MAE (Table 5).This means that the PTIM collapses into a degree-day model, leaving only the temperature-dependent term to explain melt.The temperature factor does not differentiate between snow and ice and therefore produces a poor estimate of ablation (Fig. 3a).The same result is obtained when the PTIM is tuned to ablation stake data from North and South Glaciers.
The behavior of the PTIM for Donjek Range conditions is the opposite of that documented by Pellicciotti et al. (2008) for Norte Glacier in the semi-arid subtropics of Chile.When the PTIM was applied to Norte Glacier, the temperature melt factor (T F ) became zero or negative and the model was thus entirely dependent on shortwave radiation through S RF (Pellicciotti et al., 2008).To explore whether the collapse of the PTIM to a degree-day model is an artifact of our tuning method using cumulative ablation, we instead tuned the PTIM with point-scale hourly melt rates produced by the IEBM.This was accomplished by maximizing E, with timestep j replacing location i in Eq. ( 18).This tuning method for the PTIM is similar to that of Pellicciotti et al. (2005).The best fit model using this second tuning method also produced a poor estimate of distributed ablation (Fig. 3a).
To investigate whether our implementation of the PTIM was somehow flawed, the IEBM was run with inputs appropriate to Haut Glacier d'Arolla in the study of Pellicciotti et al. (2005) (Table 4).The two methods of tuning the PTIM were then repeated with these inputs.Figure 3 shows that the PTIM produces much better estimates of ablation for Haut Glacier d'Arolla conditions than for Donjek Range conditions.The success of the PTIM in the former case is consistent with the findings of Pellicciotti et al. (2005) and Carenzo et al. (2009).
Sensitivity tests were conducted to attempt to better understand the difference in PTIM performance for these two applications.No single parameter can explain the tendency for S RF to become zero under Donjek Range conditions.However, decreases in snow albedo below 0.7 or increases in cloud fraction (above 0.4) or atmospheric refection (above 0.6) will permit non-zero values of S RF in conditions otherwise resembling those of the Donjek Range.Southerly aspects also produce a small departure of S RF from zero.Further, no single parameter alteration from the Haut Glacier d'Arolla conditions in Table 4 produces S RF = 0. Based on the results above for Donjek Range conditions, we have excluded the PTIM from the remainder of our comparative analysis.

Parameter sensitivity
The results of the parameter sensitivity tests are thematically grouped into three figures.Parameter sensitivity to glacier geometry (aspect, slope, and minimum and mean elevations) is illustrated in Fig. 4, to snow and ice albedo in Fig. 5, and www.the-cryosphere.net/5/1011/2011/The Cryosphere, 5, 1011-1028, 2011 to meteorological quantities (mean annual temperature, temperature lapse rate and wind speed) in Fig. 6.Two of the parameter sensitivity tests conducted are not shown (maximum glacier-surface elevation and initial snowpack depth) due to the small sensitivity of melt-model parameters to these quantities.
Glacier geometry affects the surface energy balance in a number of ways.Changing the aspect and slope of the glacier modifies the intensity of shortwave radiation reaching the glacier surface.Changing the minimum elevation of the glacier changes the proportion of the glacier subject to warmer temperatures according to the prescribed temperature lapse rate.Modifying the mean elevation of the glacier changes the magnitude of the turbulent heat fluxes as the glacier is moved between a warmer environment with a thicker atmosphere and a colder environment with a thinner atmosphere.One can see from Fig. 4 that DDF snow (CTIM) is strongly affected by aspect and weakly effected by the other geometric attributes.DDF ice is most sensitive to the mean and minimum glacier surface elevations.Of the HTIM parameters, r snow is most sensitive to glacier aspect (Fig. 4b), while r ice is most sensitive to mean surface elevation (Fig. 4k).The equifinality of the HTIM model is evident in Fig. 4; MF switches from being zero to non-zero under certain conditions, affecting the magnitude of the other two melt-model parameters.The parameters for the SEBM show discernible sensitivity to glacier aspect (Fig. 4c) and mean surface elevation (Fig. 4l) but of a much lower magnitude than the other models.Glacier slope (Fig. 4d-f) has little effect on the magnitude of the melt-model parameters in these simulations where the default aspect is east.The effect of slope is more pronounced for simulations conducted with a southerly aspect.For southerly aspects, the cosine of the solar azimuth (corrected for aspect) is close to unity when the sun is high in the sky.Aspect is expected to have diminishing control on melt model parameters at latitudes higher than that of our study site.During summer at high latitudes, in addition to the long hours of daylight, there is only a small difference between daily maximum and minimum solar zenith angle.
Surface albedo controls the net shortwave radiation received by the glacier and therefore has a straightforward effect on the energy balance.Figure 5 shows the sensitivity of DDF snow and DDF ice (CTIM) to the respective values of snow and ice albedo (Fig. 5a, d).DDF ice goes to zero for snow albedo values where melt is not sufficiently intense to raise the snow-line above the minimum elevation of the glacier (Fig. 5a).The strong response of CTIM parameter values to albedo has been previously shown by Arendt and Sharp (1999).In the HTIM simulations, MF is only significant for snow albedo less than 0.675 (Fig. 5b).Parameter r ice declines smoothly with increasing ice albedo (Fig. 5e), but is only non-zero for an intermediate range of snow albedo values (Fig. 5b).Parameter r snow exhibits a non-monotonic response to snow albedo that may be a function of the equifinal-ity of the HTIM.In the SEBM, C 0 is sensitive to snow albedo (Fig. 5c) while neither C 0 nor C 1 respond to changes in ice albedo (Fig. 5f).High surface albedo produces long periods where no melt occurs, this in turn allows the surface temperature the drop below zero reducing outgoing longwave radiation.It is likely that C 0 is responding this reduction in outgoing longwave radiation.
Meteorological conditions control the input of energy to the glacier system and therefore strongly influence the energy balance.Temperature appears in each of the energy balance terms, except incoming shortwave radiation, establishing a physical basis for the correlation between temperature and ablation exploited by temperature-index models (Ohmura, 2001).CTIM degree-day factors have themselves been shown to be a function of temperature by Braithwaite (1995).Wind speed affects the energy balance by enhancing or diminishing the turbulent energy fluxes.Figure 6 shows that DDF ice in the CTIM is a strong function of mean annual temperature (Fig. 6a) and the treatment of temperature lapse rates (Fig. 6d, g).DDF ice varies strongly with the true lapse rate, whether or not the true lapse rate is known to the model.In contrast, DDF snow exhibits little sensitivity to mean annual temperature or temperature lapse rates.Both CTIM parameters are sensitive to wind speed above different thresholds (Fig. 6j).The HTIM radiation factors r snow and r ice respond similarly to DDF snow and DDF ice , respectively, in this series of tests.The HTIM melt factor, MF, shows little sensitivity to the variables tested here.The SEBM parameters, C 0 and C 1 , have a discernible response to wind speed (Fig. 6l) and lapse rate (Fig. 6f, i).
The results of the lapse rate tests present an interesting dilemma.Parameters in the CTIM and HTIM models (particularly DDF ice and r ice ) vary with true lapse rate, whether or not the melt models themselves employ this lapse rate (Fig. 6g, h, tests labelled "variable") or use a standard lapse rate (Fig. 6d, e, tests labelled "fixed").In the latter case, differences between the true and prescribed (standard) lapse rates are absorbed into the tuned parameter values.This result hints that using an assumed lapse rate would not reduce the transferability of these models for the environmental conditions represented in these tests.
Overall, the sensitivity tests above suggest that the CTIM and HTIM might exhibit poor transferability where there are large variations in glacier aspect, mean glacier elevation, albedo, mean annual temperature, temperature lapse rates, and wind speed.The SEBM, as implemented, should be more transferable unless there are large variations in snow albedo or wind speed.

Parameter transferability
Table 3 reports the control run parameter values for each glacier, year and melt model, while Table 6 shows RMSE values for the control runs relative to the ablation stake measurements for each of the melt models (plus the REBM).For  In 9 of 16 temporal transferability tests (Table 6) the results more closely resemble the control runs than do those of any of the other transferability tests.This is particularly true for the REBM where all of the temporal transferability tests are close to the control run, though untrue for the SEBM where none of the temporal transferability tests are closest to the control.The results are highly variable for the spatial and spatial-temporal transfer tests but these transfers frequently produce much larger errors than the control runs.The comparison between model results in Table 6 demonstrates, within our limited data-set and bearing in mind the different tuning methods between the REBM and other models, that the REBM produces more consistent results than the temperature-index models and the SEBM.The SEBM performs inconsistently in these tests with real data, in some cases producing results better than the REBM but in other cases producing results poorer than the temperature-index models.The HTIM yields slightly more consistent results than the CTIM, but occasionally performs more poorly than the CTIM (see North Glacier temporal transferability tests).
An assessment of model performance for the master control runs is included in Table 6.In 15 of 16 cases, these simulations unsurprisingly produce higher RMSE values (from 0.01-0.29 m) than the control runs using local parameters.RMSE values for the master control runs tend to be closer to those for the local control runs for the energy balance models (SEBM and REBM), though large variations are evident.In 12 of the 15 cases above, the RMSE values for the master control runs lie between those of the local control runs and those of the spatial transferability tests.This result is intuitive considering that the parameter values derived in the master runs implicitly contain information about both glaciers.The REBM produces the best results within the master control runs for the 2008 tests and the SEBM produces the best results for the 2009 tests.
Parameters in the SEBM exhibited the least sensitivity to variations in hypothetical environmental conditions, yet do not consistently exhibit greater transferability than the parameters of the temperature-index models.This result may be related to the synthetic nature of the sensitivity tests where the SEBM knows exactly the incoming shortwave radiation and albedo, while these quantities are associated with uncertainty (particularly the evolution of snow albedo) (Gardner and Sharp, 2010) in the transferability tests with real data.

Conclusions
For the foreseeable future temperature-index models will continue to be used for melt modelling applications, due to their low data requirements and often acceptable model performance.This use should be tempered by the knowledge that for certain environments these models may exhibit poor transferability.That is, using melt-model parameters derived for one glacier at other sites can lead to large errors in estimated ablation.Driving an energy balance model with idealized inputs and tuning the melt models of interest to the output can be a useful means of exploring the sensitivity of melt-model parameters to variations in glacier geometry, surface characteristics and meteorological conditions.This approach may provide insight into the transferability of various model parameters to different environments, but has limitations as demonstrated by the performance of the SEBM with real and synthetic data.
In this study melt-model parameters showed only small differences when tuned to cumulative ablation using E, RMSE or MAE, in contrast to using a weighted R 2 .Such small differences will not significantly affect the sensitivity or transferability of the models.When the temperatureindex model of Pellicciotti et al. (2005) is tuned to output of the EBM forced with either idealized inputs (for synthetic glaciers with northerly, easterly or westerly aspects) or real ablation data from our study area in the St. Elias Mountains, the solar radiation factor S RF becomes zero and the model collapses into a classical temperature-index model.We hypothesize that some combination of study-area characteristics and meteorological conditions produces this model behaviour, but further work is required to more precisely assess the importance of each of these attributes.
In tests where the melt models were tuned to output from the idealized EBM, the temperature-index model parameters were generally sensitive to glacier aspect, mean surface elevation, albedo, mean annual temperature, temperature lapse rate and wind speed.In transferability tests using real data from two glaciers and two melt seasons, the temperatureindex models (HTIM and CTIM) produced errors up to eight times larger than their respective control runs; the simplified energy balance model (SEBM) produced errors up to six times larger than its control runs in the same tests.The simplified energy balance model more often than not produced the best parameter transferability, but its poorest transfers are just a poor as those from the temperature-index models.When compared to parallel transferability experiments reported by MacDougall and Flowers (2011) the full energybalance model produces better transfers nine times out of twelve.Using "master" parameter values derived from both glaciers and both years usually produced errors in simulated ablation higher than those obtained with locally derived values (local control runs) but lower than those with parameters tuned for the other site.The results of these transferability tests contrast with those from previous studies, where temperature-index models were found to be transferable under most conditions explored (Shea et al., 2009;Carenzo et al., 2009).This contrast could be an artifact of the small sample size in this study (two glaciers and two melt seasons) or other factors yet to be elucidated.We suggest that caution should be observed when extending the use of melt models beyond the locations where they were developed, particularly if the data are limited in spatial or temporal scope.

Fig. 1 .
Fig. 1.Study area.(a) Donjek Range between Kluane and Kaskawulsh Glaciers.(b) Surface contour map of South Glacier with locations of ablation stakes, temperature microloggers, and AWS.(c) As for (b) but for North Glacier.

Fig. 2 .
Fig. 2. Diagram of transferability tests.Arrows indicate the possible spatial and temporal transfer of parameter values between data sets.

Fig. 3 .
Fig. 3. Ablation as a function of elevation computed with the temperature-index model of Pellicciotti et al. (2005) tuned to both cumulative ablation and hourly ablation rates as simulated by the IEBM.The energy balance model forced with idealized inputs (IEBM) is also shown.(a) Inputs intended to represent Donjek Range glaciers (DRG in Table 4).(b) Inputs intended to represent Haut Glacier d'Arolla (HGA in Table4).

Table 1 .
Instrumentation deployed at AWS locations on North and South Glaciers.Manufacturer documentation is the source of instrument precision.Precipitation for South Glacier was measured 500 m from AWS.

Table 2 .
Mean values of meteorological variables and ice albedo for North and South Glaciers as measured at AWS locations from 19 May to 28 July in 2008 and 2009.Winter balances are also shown.S08 is South Glacier 2008, S09 is South Glacier 2009, N08 is North Glacier 2008 and N09 is North Glacier 2009.

Table 3 .
Parameters used for or derived in the control runs where melt models are tuned to real data for South Glacier (S) and North Glacier (N) for 2008 (08) and 2009 (09).In the master control run (M) data from both glaciers in both years are used together to derive model parameters.
MacDougall, 2010 andacier melt model parameter sensitivity and transferability to cumulative melt data, the REBM is tuned to independent measurements of snow albedo and aerodynamic roughness length (seeMacDougall, 2010 and MacDougall and Flowers,  2011for details).

Table 5 .
Best fit parameters for CTIM, HTIM, PTIM and SEBM found by optimizing slope weighted coefficient of determination (wR 2 ), Nash-Sutcliffe efficiency criterion (E), mean absolute error (MAE), and root mean square error (RMSE) with respect to cumulative melt simulated with the IEBM for synthetic glaciers oriented in the four cardinal directions.
Melt-model parameter sensitivity to glacier aspect (first row), surface slope (second row) and minimum (third row) and mean (fourth row) elevation.Minimum elevation and glacier length are manipulated together in order to preserve glacier slope.For the mean elevation test, initial snow depth is held constant.Note the tangential scale for slope.bothglaciers, the model with the lowest RMSE is the HTIM for the 2008 simulations and the SEBM for the 2009 simulations.The highest RMSE is produced by the SEBM for both glaciers in 2008 and the CTIM for North Glacier 2009.