Reconstructing glacier-based climates of LGM Europe and Russia - Part 1: Numerical modelling and validation methods

. The mountain environments of mid-latitude Europe and Arctic Russia contain widespread evidence of Late-Quaternary glaciers that have been attributed to the Last Glacial Maximum (LGM). This glacial-geological record has yet to be used to quantitatively reconstruct the LGM climate of these regions. Here we describe a simple glacier-climate model that can be used to derive regional temperature and precipitation information from a known glacier distribution. The model was tested against the present day distribution of glaciers in Europe. The model is capable of adequately predicting the spatial distribution, snowline and equilibrium line altitude climate of glaciers in the Alps, Scandinavia, Caucasus and Pyrenees Mountains. This veriﬁcation demonstrated that the model can be used to investigate former climates such as the LGM. Reconstructions of LGM climates from proxy evidence are an important method of assessing retrospective general circulation model (GCM) simulations. LGM palaeoclimate reconstructions from glacial-geological evidence would be of particular beneﬁt to investigations in Europe and Russia, where to date only fossil pollen data have been used to assess continental-scale GCM simulations.


Introduction
To provide confidence in climate predictions made using general circulation models (GCMs) it is important to compare their predictions of past climates with records of past climates. GCMs require observations and measurements for model inputs and boundary conditions as well as information against which the model can be tested. The Last Glacial Correspondence to: M. J. Siegert (m.j.siegert@ed.ac.uk) Maximum (LGM) is the most recent prelonged cold phase in the Earth's history (e.g. Siegert 2001). Owing to the different nature of the climate and relative abundance of preserved evidence for climate change the LGM is a popular time period for testing the ability of GCMs to simulate past climates (e.g. the Palaeoclimate Modelling Intercomparison Project (PMIP) (Joussame and Taylor, 1995) and PMIP2 collaborative projects (Harrison et al., 2002). To date the only continental-scale proxy LGM climate reconstructions used to assess GCM simulations of Europe and Russia have been derived from fossil pollen data (Peyron et al., 1998;Tarasov et al., 1999;Kageyama et al., 2001;Jost et al., 2005;Wu et al., 2007;Ramstein et al., 2007). It is important to try and use a multi-proxy approach, such as established in the tropics (e.g. Farrera et al., 1999) when assessing GCM model output for the following three reasons. First, a single proxy source may not provide a complete climate reconstruction. Individual proxy records will primarily reflect the aspects of the climate to which they are most sensitive; plants (and therefore pollen) will most reliably reflect "bioclimatic" variables (e.g. temperature of coldest month, or seasonal distribution of precipitation) (Prentice et al., 1992), rather than "traditional" climate variables (e.g. mean annual temperature or annual precipitation). Second, methodological limitations may create errors in the reconstructed climate signal. The method used by Peyron et al. (1998) and Tarasov et al. (1999) assumed that the change in vegetation distribution between the present day and LGM reflected a change in climate alone. Modelling studies (e.g. Jolly and Haxeltine, 1997;Harrison and Prentice, 2003) and laboratory studies (Cowling and Sykes, 1999) have shown that the distribution of LGM vegetation is affected by the reduced atmospheric CO 2 concentration during the LGM. The omission of this factor from the Peyron et al. (1998) and Tarasov et al. (1999) reconstructions means that the LGM precipitation anomaly is over-estimated (Cowling and Sykes, 1999). Third, a multi-proxy approach allows regional trends reconstructed within a single proxy to be corroborated (e.g. Farrera et al., 1999). This is important because the coarse resolution of GCMs prevents them from simulating local scale factors that influence the climate signal recorded at individual proxy sites.
Consequently, there is a need for new continental-scale proxy LGM climate reconstructions across Europe and Russia which can contribute to the calibration of present and future GCMs. Glaciers can be used as indicators of environmental change; the spatial distribution of glaciers is, to a first order, a function of precipitation and temperature conditions. The climate conditions required to maintain individual glacier mass balance have been modelled in a variety of ways (e.g. Oerlemans, 1991;Hock, 1999;Braithwaite and Zhang, 2000;Bassford et al., 2006). In order to compare the results of large-scale climate models with those derived from glaciers, a method is needed which can extract regional information concerning precipitation and temperature from mass balance models.
This paper is the first of a trilogy aiming to demonstrate how the extent of former glaciers and ice caps can assist large-scale palaeoclimate reconstruction. In this first paper, a modelling approach designed to characterise the regionalscale relationship between climate and glaciated regions is presented. The model is tested by applying it to the present day, to see whether it can predict the known extent of glaciers given modern accumulation and temperature records. The result is a model capable of determining the climate required under a given distribution of glaciers (e.g. at the LGM and in the future). In the second paper (Allen et al. 2007a) the model is used to understand the relationship between temperature and precipitation across Europe at the LGM, given the published extent of glaciers at this time (Ehlers, J. and Gibbard, 2004). In the third paper (Allen et al., 2007b), this relationship is tested against reconstructions of LGM climate based on GCM modelling and on pollen assemblages. The trilogy provides a complete methodology for how past glaciers can be used in other parts of the world, and in future as both modelling and our knowledge of past glacier extents improves.

The degree day model
The mass budget and extent of glaciers are determined by the climate and the characteristics of ice (see Paterson, 1994 for a full review). The geological record demonstrates that glaciers are sensitive to changes in climate (e.g. Ehlers and Gibbard, 2004). It is on this premise that glacial-geological evidence has been widely used to make inferences about past climates (e.g. Leonard, 1989;Kull et al., 2003;Mark et al., 2005). A degree day model (DDM) was used to calculate ablation at the glacier surface in this study. This approach uses the sum of positive air temperatures (T + ) to calculate melting (M) during a given time period ( t(d)), divided into n time intervals, the factor of proportionality is controlled by the degree day factor (DDF) expressed in mm d −1• C −1 (Hock, 2003): It is usual to use different degree day factors for snow and ice surfaces to reflect the lower albedo and higher ablation rates of ice compared to snow (Hock, 2003). Surface accumulation is calculated using a temperature threshold to divide precipitation (P * ) into rainfall or snowfall: The mass balance model was used to simulate mass balance over a pre-defined glacier geometry using the principles of static mass balance sensitivity. This approach assumes the glacier geometry remains fixed and does not explicitly calculate glacier flow. The response of the glacier to climate is represented by changes in the mass balance profile from the fixed glacier geometry (e.g. van de Wal and Oerlemans, 1994;Fleming et al., 1997). Static sensitivity experiments on palaeo-glaciers assume steady-state conditions. The mass balance model is tuned until cumulative surface mass balance is zero representing equilibrium in the glacier climate system (e.g. Hostetler and Clark, 2000).

Numerical details
It is not possible to derive melt factors for LGM glaciers, therefore the model was parameterised using melt factors measured over present day glaciers, and it is assumed these values adequately represent the LGM climate-glacier relationship. The average melt factors for Scandinavian and Alpine glaciers from Braithwaite and Zhang (2000) are 4.3 mm d −1• C −1 and 6.5 mm d −1• C −1 for snow and ice, respectively; these values were used in this study and the melting threshold was set at 0 • C. For mid to high latitude glaciers the precipitation threshold is usually between 0 • C and 2 • C (e.g. Bassford et al., 2006); a value of 1 • C was used to incorporate the occurrence of snowfall above 0 • C. Rainfall and meltwater were assumed to runoff the glacier surface in the model and make no contribution to net accumulation via refreezing or superimposed ice formation.
To ensure the numerical stability of the mass balance calculations each simulation was initiated with a default snow surface. Ablation and accumulation (Eqs. 1 and 2) were calculated on an hourly timestep and the model was run for one year starting on 1st September (Julian Day 244) (assumed to be the start of the winter accumulation season). This allowed the development and melting of the snowpack during the winter and spring, respectively. Once melting had started in the spring, the equivalent melt from each time step was removed from the snowpack. If the snowpack was melted away completely the model switched to melting the ice surface.

Applicability of the model
The simplicity and requirement for only two meteorological parameters (temperature and precipitation) mean that DDMs have been widely used in palaeoclimate modelling studies (e.g. Hostetler and Clark, 2000;Kull and Grosjean, 2000;MacGregor et al., 2000). The simplicity of the DDM allows an ease of application, especially in regions where data are limited; the trade off is that there are limitations to what they can achieve. Using fixed degree day factors only bulk "average" conditions can be estimated and local-scale glacial processes will not be captured (Hock, 2003). DDMs are insensitive to changes in the style of seasonality, specifically the winter season; once air temperature drops below the melting threshold ablation will cease and the magnitude of the negative temperature is not considered. Static-mass balance sensitivity analyses have limited applicability when studying climate change in the recent past because the dynamic response of mountain glaciers is short (10 to 10 2 years), therefore an appreciation of changes in hypsometry is required to fully understand the glacier response to the climate signal. For studies investigating longer-term climate variation (10 3 to 10 4 years) it can be assumed that glacier changes are a response to longer term mean climate forcing making the assumption of steady-state more plausible (Seltzer, 1994).

Basic modelling procedure
The basic modelling procedure is as follows. The model attempts to predict the extent of glacierisation, which can be compared with modern and ancient glacier coverage data. The model solution is non-unique, however, as there are a range of temperature/precipitation arrangements that lead to very similar glacier sizes. For the modern case, we can use present-day climate data to test whether the model is capable of predicting modern glacier sizes. For the ancient case, however, the model simply provides all the possible arrangements of precipitation and accumulation responsible for glacierisation (Allen et al., 2007a), which can be tested against other reconstructions of palaeoclimate (Allen et al., 2007b). For the LGM reconstructions, we use present-day climate data as our starting position. We then modify temperature and precipitation and (lapse rates) accordingly to get the LGM distribution of glaciers. In doing so, we show the possible adjustments needed to the modern climate regime to get the glaciers known to have existed in the past.

Input data
To calculate glacier cumulative mass balance the DDM requires a hypsometric profile (i.e. the spatial distribution of the glacier as a function of altitude). In this study the spatial geometry and altitudinal profiles of the present day or LGM glaciers were reconstructed separately and combined to produce the hypsometry used by the DDM. The USGS "gtopo30 arcsec" DEM (USGS, 1996) was used to provide the altitudinal component of the present day and LGM climate reconstructions. The resolution of this DEM (∼1 km) provides a good representation of the broad scale relief and altitude range within the upland regions glacerised now and glaciated at the LGM.
Details of the spatial geometry used to represent the present day glaciers and LGM glaciers are described in Sect. 3.2 of this paper and in Allen et al. (2007a), respectively.
The high resolution (10 latitude/longitude) monthly CRU2.0 climate dataset, created by the Climate Research Unit (CRU), University of East Anglia, was used to represent the present day climate baseline from which LGM climate anomalies would be derived. This data set was constructed using a thin-plate spline interpolation for the period 1961-1990. The spline interpolation is a three-dimensional (i.e. altitude sensitive) interpolation (Hutchinson, 1999). New et al. (2002) provide a full description of the CRU2.0 climate dataset, which has three advantages of relevance to this study. First, the dataset enables all the simulations (including those in Allen et al., 2007a, b) to be driven with meteorological data from the same source constructed using a consistent methodology. Second, the dataset represents a 30-year average climate; a single-year climate record may not necessarily be representative of a mean present day climate. Third, the individual meteorological variables are accompanied by an uncertainty (New et al., 2002) enabling the sensitivity and reliability of model results to be tested against the uncertainty of the input data.
The present day climate used to drive the model was based on the mean monthly temperature ( • C), mean monthly diurnal temperature range ( • C) and monthly precipitation totals (mm mo −1 ) from the CRU2.0 climate dataset (New et al., 2002). These variables are presented in the dataset on a monthly resolution; they were downscaled to the diurnal climate required by the DDM as follows: it was assumed that precipitation rates were constant throughout each month and hourly precipitation (P h ) was calculated from the CRU2.0 monthly precipitation total as (P CRU ) and days in the month (dpm): Alps 5327   The hourly temperature (T h ) was calculated from the mean monthly temperature (T mo ) and diurnal temperature range (T r ) using a cosine function similar to that used by Bassford et al. (2006): Within each month the mass balance totals simulated over the diurnal cycle where scaled up to form the monthly mass balance total. The CRU2.0 climate dataset was downscaled to each cell in the USGS DEM (USGS, 1996) using temperature and precipitation lapse rates. Owing to the absence of field measurements that could be used to prescribe site specific lapse rates they were treated as unknowns in the modelling experiments. To encompass all possibilities, a suite of 189 lapse rates was used to represent temperature lapse rates ranging from 0 • C/km to 10 • C/km, and precipitation lapse rates ranging from 0 mm/100 m to 80 mm/100 m (this range is similar to the range of published precipitation lapse rates across Europe, e.g. Sevruk, 1997). To be clear, we use only one lapse rate per experiment; we do not apply spatially varying lapse rates, though such application might be of interest in future work. The precipitation lapse rate was used to adjust the annual precipitation total and the resulting change in precipitation applied evenly across the year. It is acknowledged that the downscaling of the CRU2.0 climate is an extrapolation of the dataset and will create model climates at the DEM resolution that were not used in the creation of the CRU2.0 dataset. Furthermore it is accepted that the lapse rates are be-ing used in a purely pragmatic modelling context and are not attempting to simulate the physical processes by which lapse rates occur.

Model test data
Before being applied to retrospective climate reconstructions (Allen et al. 2007a, b), the ability of the DDM to characterise the regional scale glacier-climate signal of five currently glacierised regions in Europe, which are characterised by a small glacierised extent but a high number of individual discrete valley and mountain glaciers, is tested (Table 1). The only dataset containing the level of detail to adequately describe these regions is the World Glacier Inventory (WGI) (National Snow and Ice Data Centre, 1999). To make the ASCII formatted WGI data compatible with the DDM results it was necessary to convert it into a grid format. Whilst the WGI data describes the size, altitude range, and total area of individual glaciers, it provides no hypsometric data. As a result it is impractical to construct individual glacier profiles at the DEM resolution, especially for glaciers greater than 1 km 2 . The WGI data were converted to a grid with the same resolution as the CRU2.0 climate dataset. Each glacier was attributed to a grid cell using the latitude and longitude attributes. For each cell the contributing glacier data were used to construct a total glacierised area (Fig. 1) and average snowline, maximum, minimum and mean altitude. It is noted that the WGI descriptive data (snowline altitude, maximum, minimum and mean altitude, and glacier area) for the glacierised regions of Europe are incomplete (Table 2). Therefore the cell characteristics derived from the combined WGI data may not wholly reflect the glacial characteristics of each ∼20 km cell.
The model was used to simulate the annual mass balance for all the DEM cell within the model domains defined for each region (Table 3) which were scaled up to the resolution of the CRU 2.0 dataset for comparison. CRU cells containing DEM cells with positive annual mass balance were assumed to be glacierised, and conversely cells containing only negative annual mass balance DEM cells were assumed to be non-glacierised. Four types of result could be predicted by the DDM when compared to the WGI dataset. Type one, correctly predicting the location of a WGI glacierised cell. Type two, correctly predicting the location of a WGI nonglacierised cell. Type three, predicting a WGI glacierised cell to be non-glacierised. Type four, predicting a WGI non-glacierised cell to be glacierised. A simple cost function was used to optimise the lapse rate combination that minimised the difference between the model predictions and WGI dataset. The cost function (CF) calculated the number of type one (A) and type two (B) results and compared them to the number of glaciated (A ) and non-glaciated cells (B ) in the WGI dataset, Eq. (5). The cost function returns a value between 0 and 1, with one indicating a perfect prediction of the WGI dataset by the model.  Kotlyakov and Krenke (1982) and Leonard (1989). An alternative dataset of measured ELA climates (Ohmura et al., 1992) is plotted as a compariosn. Apart from two Alpine ELA climates measured by Ohmura et al. (1992) the agreement between the two independently derived datasets is good and provides confidence in the use of the ELA climate "envelope" as a method for assessing DDM predictions of ELA climate.
A comparison of the summer temperatures and winter precipitation at the ELA (Fig. 2) were used to assess the glaciological and climate conditions simulated by the model. The model ELA was calculated as a function of the altitude and mass balance between neighbouring DEM cells: where, E 1 and E 2 are the elevation of neighbouring DEM cells with positive and negative mass balance, respectively,   and bn 1 and bn 2 are the annual mass balance of the positive and negative DEM cells, respectively (Oerlemans, 1991). This comparison assumed that the ELA and snowline on the modelled glaciers are at the same altitude. This is reasonable because they are generally found at similar altitudes on temperate mountain glaciers (Benn and Evans, 1998), although it is acknowledged that they are different glaciological parameters. The envelope of present day ELA climates is based on data from 32 glaciers (Kotlyakov and Krenke, 1982;Leonard, 1989).

Experiment one: spatial distribution of glacierised and non-glacierised regions
The aim of this experiment was to quantify the ability of the DDM to simulate the known distribution of present day glaciers in the five model regions. Across the suite of climate lapse rates the cost function results are more sensitive to the temperature lapse rate than the precipitation lapse rate. Using small temperature lapse rates the DDM simulates negative mass balance in all cells, as temperature lapse rate increases the number of correctly predicted glacierised cells (and cost function) increases (Fig. 3). At very high temperature lapse rates the percentage of correctly predicted non-glacierised cells (and cost function) starts to decline. In Southern Scandinavia and Caucasus Mountains there are multiple cost function optima, with the same cost function value, but different spatial predictions of glacierised and nonglacierised cells (Table 4). Under optimum lapse rates the DDM predictions of regional glacierization follow the same broad pattern in all regions; predictions of non-glacierised zones exceed 90%, and predictions of glacierised cells exceed 50% (Table 4 and Fig. 4). The accumulation area ratio (AAR) of a glacier describes the proportion of an accumulation zone relative to the total glacier area. Published AAR values for mid-high latitude glaciers range from 0.5 to 0.8 (Benn and Lehmkuhl, 2000), with 0.67 being a commonly used value (Benn and Evans, 1998). The total glacial extent within each glacierised cell was estimated using an AAR of 0.67; assuming that DEM cells with positive mass balance represented the accumulation zone of the glaciers. In the five modelled regions, the within-cell glacial coverage described by the WGI is less than 5% in the majority of glacierised cells. The DDM predicts a similar extent of glacial coverage (Table 5) (Allen, 2006). The aim of this experiment was to assess if the ELA climates simulated by the DDM were compatible with the measured ELA climates described in Fig. 2. In each region the ELA climate was derived from the simulations using the optimised lapse rate combination determined in Experiment One. The model is able to predict well the summer temperatures and winter accumulation rates at the ELA in most regions (Fig. 5). One exception is in the Caucasus Mountains, where the simulated ELA climates straddle the left hand boundary of the measured ELA climate envelope. The modelled ELA climates lying outside of the ELA climate envelope are, however, no more extreme than the outliers in the Ohmura et al. (1992) dataset.

Experiment three: DDM simulated ELA estimates
The aim of this experiment was to assess the ability of the DDM to replicate the altitudinal profile of the glaciated regions described by the WGI. For each region only the cells with the most complete WGI dataset were selected for this experiment. In this experiment the optimum lapse rate combination which minimised the difference between the model ELA and WGI snowline data was used as the optimum result. In the Alps the DDM could simulate the ELA to within 100 m of the mean WGI snowline in 11 of the 12 assessed cells (Fig. 6). In other words, information on the presence of glaciers (at the regional scale) and local information on climate (CRU) is sufficient for the model to predict the ELA. In the Caucasus and Scandinavian regions the DDM ELA estimates were lower than the maximum glaciated altitude, but systematically higher than the mean WGI snowline. The discrepancy between the lowest DDM ELA prediction and the within cell mean WGI snowline ranged from 162 m to 309 m

Experiment four: sensitivity analysis
The aim of the sensitivity analysis was to investigate the extent to which the DDM predictions of present day European glaciers changed in response to first, uncertainty in the input data and second, the range of potential DDM parameterisations. Such investigation is required to fully understand how representative the climate created by the modelling approach is of the present day. The uncertainty in the CRU2.0 climate data, vertical error in the USGS DEM and eight different DDM parameterisations were combined to create a suite of sixteen sensitivity experiments ( Table 6). The same methodologies used in the first three experiments were used to determine the optimum model predictions of glacier distributions, ELA and ELA climates for each model sensitivity simulation.
In the Alps, Pyrenees, Caucasus Mountains and Southern Scandinavia the optimum model result from each model sensitivity experiment predict a distribution of glacierised and non-glacierised cells comparable to the model control runs (Fig. 7). The re-optimized lapse rate combinations do vary compared to the control experiment, and reflect that, with other things being equal, the changes to the input data or model parameters will either increase or decrease the annual mass balance simulated by the DDM. In simulations which increase the annual mass balance (TEMP-, RANGE-, PPT+, ELEV+, DDF-1, DDF-2, THOLD-2, and THOLD-4) the re-optimised temperature lapse rate is decreased compared to the control. The inverse occurs in the remaining experiments which decrease the annual mass balance total. Northern Scandinavia is the most sensitive region to the uncertainty in input data and potential range of parameterisations; this is reflected in changes to the correct prediction of glacierised cells of up to 18% from the control. In all regions the sensitivity experiments cause only small scale changes in the DDM predictions of within-cell glacier coverage, ELA estimates and ELA climates. These variations are not significant enough to change the regional trends present in the control experiments.

Discussion
A key characteristic of glacierised regions is the distribution of glaciers and surrounding non-glacierised zones. Using the optimum lapse rate combinations the model predicted >90% of the non-glacierised cells and >50% of glacierised cells in all modelled regions. The structure of the cost function results from experiment one and sensitivity analyses indicate that the model predictions were statistically the best achievable results for all regions except Northern Scandinavia. At high lapse rates, the presence of type four results (WGI nonglacierised cells predicted by the DDM as glacierised) show that if the lapse rate domain had included larger lapse rate values than those used the overall cost function would not have increased. Whilst using larger lapse rates would increase the percentage of correctly predicted glacierised cells this positive effect on the cost function would have been negated by the increasing presence of type four results. In the sensitivity analyses the different experiments changed the optimum lapse rate combination, however the predictions of glacierised and non-glacierised cells were not changed significantly from the control simulation. In Northern Scandinavia, the cost function value increases significantly in sensitivity analyses where the DDM is able to simulate more positive annual mass balance compared to the control simulation. This suggests that Northern Scandinavia has a distinct climate regime compared with other parts of Europe, which makes sense given its maritime environment.
The inability of the DDM to correctly predict higher percentages of glacierised cells is most likely to be related to the characteristics of the DDM, USGS DEM and WGI datasets. Despite the glacierised regions considered in this study containing numerous glaciers, the individual glaciers are relatively small. As such the majority of glaciers are likely to be influenced by significant local topographic or climatic factors, e.g. steep sided valleys reducing direct insolation, topographically induced precipitation, or wind blown snow. These local scale processes cannot be reproduced by the CRU2.0 and USGS DEM datasets. It is possible that some of the glaciers detailed in the WGI dataset are sustained by these processes in regions where the regional climate does not alone sustain glacierization. Many of the glaciers have a surface area that is beneath the resolution of the DEM, e.g. in the Pyrenees the largest glacier is ∼1 km 2 . In such cases the DDM will return a single mass balance value to represent the whole glacier. If this is negative the cost function would assume that the region is non-glacierised. Higher resolution (<100 m) DEMs (e.g. Shuttle Radar Topograhic Mission) are now becoming available, and would provide a more detailed model representation of the topography in mountainous regions. The application of such DEMs is currently limited in the regional scale modelling discussed in this paper owing to the resolution of available climate data, which does not contain sufficient detail to be reliably downscaled onto a DEM with a resolution <100 m. The dates of the WGI observa- Fig. 5. Climate at the ELA predicted by the DDM in the Alps using the optimum lapse rate combinations derived in Experiment One compared to ELA climates measured on present day mid-latitude mountain glaciers (Kotlyakov and Krenke, 1982;Leonard, 1989). tions, used to characterise the glaciers within the modelled regions, range from 1952 to 1983. Therefore some of the WGI data predate the CRU2.0 dataset, which represents the climatological normal 1961-1990. Global temperatures have shown a warming trend and the mass balance of European glaciers has been generally negative during the 20th Century (IPCC, 2001). It is possible that some of the smallest glaciers contained in the WGI dataset ceased to exist between 1961 and 1990. Consequently the WGI maps used to assess the DDM predictions have to be viewed as a maximum glacial characterisation of the period 1961 to 1990.
Within each model region the DDM was able to simulate a style of ELA climate that is compatible with measured ELA climates, demonstrating that the modelling approach could consistently create plausible climatic conditions over glacier surfaces. The positively skewed distribution of "within-cell" glacier coverage predicted by the DDM in all regions is in broad agreement with the style of glacier coverage described in the WGI; however there are important differences between the model results and the WGI data that must be discussed further.
In all regions except the Alps the mean, and range, of predicted "within-cell" glacial coverage are smaller than the WGI dataset. This reflected in the systematic over prediction of the ELA by the DDM compared to the mean WGI 244 R. Allen et al.: Glacier-climate model for reconstructing palaeoclimates Fig. 6. WGI (National Snow and Ice Data Center, 1999) within cell altitude distributions and DDM ELA estimates. The maximum glacierised altitude in the WGI dataset is indicated by the star, the boxplot beneath the star is the altitudinal distribution of WGI snowline measurements within the cell. For identification purposes during model simulations cells in the ∼20 km resolution model domains (see Table 3) were numbered starting from the top left hand corner of the model grid and finished at the bottom right hand corner, each row was numbered left to right. Cell numbers have been included in this diagram to help the reader compare WGI data (left hand panel) and DDM results (right hand panel) from the same cell.
snowline. As part of a study of LGM glaciers in the tropics Hostetler and Clark (2000) verified their DDM by simulating modern tropical glaciers. They used the USGS DEM and climate predictions from the GENESIS (v.2.01) general circulation model. Whilst their DDM could simulate the ELA and the mass balance gradient the spatial extent of the glaciers was over predicted by 50%. They attributed this over pre-diction to first, local scale topographic features that create favourable climatic conditions required for glaciation, and second, the size of the ablation area of tropical glaciers being beneath the resolution of the DEM. Therefore, successful simulations of the altitudinal range of the glaciers required the DDM to over predict the glaciated area. If this interpretation of the effect of the USGS DEM resolution on the DDM Fig. 7. Percentage of correctly predicted glacierised and non-glacierised cells across the suite of sensitivity analyses. See Table 6 for details of individual sensitivity experiments. simulation is correct, it potentially indicates that the modelling approach tested in these verification experiments created a climate that was either too warm or dry, and caused the DDM to simulate annual mass balance values that were too small. This would explain the ELA over-predictions, underpredictions of "within-cell" glacial coverage, and optimum lapse rate combinations. The optimum temperature lapse rates are higher than the environmental lapse rate (6.5 • C/km) which can be viewed as the idealised optimum lapse rate because it is commonly measured (Barry and Chorley, 2003) and frequently used in climate modelling studies (McGuffie and Henderson-Sellers, 1997). A climate that is either too www.clim-past.net/4/235/2008/ warm or dry will optimise at a higher lapse rate to enable the model to exaggerate the altitudinal influence on the climate.
The excellent all round results in the Alps suggest that the model climate bias is spatially variable. The Alps have a long history of both climate and glacier observations collected from a dense network of observation posts. As a result, it would be expected that the altitudinal influence on climate and glacier measurements in the Alps would be well represented in the CRU2.0 and WGI datasets, respectively. This level of detail is not available for the other modelled regions (see Figs. 1-9 in New et al., 2002). As stated in Sect. 3.3 the downscaling of the CRU2.0 dataset used in this study is an extrapolation, and therefore may create erroneous model climates in mountain regions which are not fully represented in the CRU2.0 climate (Allen, 2006).

Conclusions
A simple method by which glaciers provide climatic information at a regional scale has been outlined. The model using modern climate as an input was tested against the known record of glaciers in the Alps, Pyrenees, Scandinavia and Caucasus Mountains, and was found to be capable of predicting the distribution and characteristics of these currently glacierised regions. In the five modelled regions the DDM correctly predicted over 90% of the non-glacierised cells, and between 50% and 87% of the glacierised cells (Table 4), furthermore the distribution of glacierised cells and the within cell glacial extents predicted by the DDM were in good agreement with the WGI data ( Fig. 4 and Table 5). The ELA climates predicted by the DDM correlate with ELA climates measured on European glaciers (Ohmura et al., 1992) (Figs. 2 and 5). In the Alps, where the glacier data are most reliable and the meteorological network is dense, the DDM was able to reliably simulate the average snowline