Snow glacier melt estimation in tropical Andean glaciers using artificial neural networks

Snow and glacier melt (SGM) estimation plays an important role in water resources management. Although melting process can be modelled by energy balance methods, such studies require detailed data, which is rarely available. Hence, new and simpler approaches are needed for SGM estimations. The present study aims at developing an artificial neural networks (ANN) based technique for estimating the energy available for melt (EAM) and SGM rates using available and easy to obtain data such as temperature, short-wave radiation and relative humidity. Several ANN and multiple linear regression models (MLR) were developed to represent the energy fluxes and estimate the EAM. The models were trained using measured data from the Zongo glacier located in the outer tropics and validated against measured data from the Antizana glacier located in the inner tropics. It was found that ANN models provide a better generalisation when applied to other data sets. The performance of the models was improved by including Antizana data into the training set, as it was proved to provide better results than other techniques like the use of a prior logarithmic transformation. The final model was validated against measured data from the Alpine glaciers Argenti ère and Saint-Sorlin. Then, the models were applied for the estimation of SGM at Condoriri glacier. The estimated SGM was compared with SGM estimated by an enhanced temperature method and proved to have the same behaviour considering temperature sensibility. Moreover, the ANN models have the advantage of direct application, while the temperature method requires calibration of empirical coefficients.


Introduction
Glaciers could be considered as the most important water reservoirs, since they represent about 68% of the total fresh water available (Shiklomanov and Roda, 2003). Previous studies state that snow and glacier melt (SGM) is of fundamental importance for the present and future water scenarios in snow-fed and glacier-fed basins (Kure et al., 2012a;Jansson et al., 2003), but most of them are located in the poles far from human activities; only mountainous glaciers are located in human populated continental areas. Mountainous glaciers could be considered the world's virtual water towers assuring year round water flow for the main rivers, and its melting may lead to water shortage for millions of people. Unfortunately, most of the mountain glaciers are melting quite rapidly, a fact that may lead to serious social tensions related to water. Hence, it is important to understand glacier dynamics in order to analyse possible future water scenarios.
There have been many studies about glaciers and snowfall, both at global and local scales. Radic and Hock (2010) applied a statistical method to estimate global glacier volume and states that corresponds to a sea level equivalent of 0.7 m. Avian and Bauer (2006) monitored Pasterze glacier with laser scanning technique and detected three zones of collapsing ice body. Huss et al. (2010) analysed the spatial distribution of Switzerland glaciers by relating glacier surface elevation change as a response to mass balance change. Koboltsching and Schöner (2011) investigated the contribution of glacier melt to total river run off in the Austrian Alps. Also different measures to prevent melting like covers, water injection or snow compaction were tested at field locations (Olefs and Fischer, 2008). Nevertheless, the above mentioned studies were applied to high latitude locations where climatic conditions are different than tropical latitude places like the Andes.
The dynamics of the tropical Andes is quite different than the Alpine glaciers; while the Alps experience a long accumulation period in winter, the Andes experience permanent ablation throughout the year (Coudrain et al., 2005). In both inner and outer tropics the daily variation is greater than the seasonal one, which is typical behaviour of tropical latitudes (Mote and Kaser, 2007). The tropical Andean glaciers used to cover over 2940 km 2 , but suffered a strong retreat that decreased its area to 2493 km 2 by 2002 and caused some small glaciers such as Chacaltaya or Cotacahi to disappear (Vergara et al., 2007), with serious consequences. For instance, the area around Cotocachi not only experienced a decrease in agriculture and tourism activities but also more and worse water conflicts are expected over time. Important Bolivian cities such as La Paz and Cochabamba already faced serious social tensions categorised as emblematic in global water debates (Laurie and Crespo, 2007). Thus, a better water resources management is an important goal. SGM estimation is an important and necessary tool for achieving such goal.
Analytical methods for estimating SGM may be divided into temperature index models and energy based models. Although temperature index models are a simplification of complex processes that would be better described by energy balance, many studies found a high correlation between melt and air temperatures . Temperature models relate the amount of melting to a degree day factor and to either the sum of the positive temperature or the mean daily temperature. Sometimes temperature index models use a base temperature that might be below the freezing temperature (Debele et al., 2009). Since they have the advantage that temperature is an easy to measure data, they are popular and used in many studies (Kure et al., 2012b;Hock, 1999;Jost et al., 2012;Biggs and Whitaker, 2012). Nowadays, some hydrological models include the option of snow/ice melting by using temperature based equations that were used in different studies (Abbott et al., 1986;Scharffenberg and Fleming, 2010;Wang et al., 2010;Tahir et al., 2011;Bocchiola et al., 2011). Hirabayashi et al. (2010) estimated global glacier mass balance using the global glacier model HYOGA that uses a day degree approach; nevertheless, such model was not able to simulate the temporal variation in South American glaciers. Pelliccioti et al. (2005) developed an enhanced temperature index model that combines temperature with shortwave radiation for estimating SGM in Switzerland (Carenzo et al., 2009).
Earlier studies showed that just air temperature is not enough for predicting snowmelt (Zuzel and Cox, 1975). It is important to consider that temperature of matter is just a property that represents the relation between the heat added to a body and its change in entropy (King, 2005), thus the external heat added to a given body, in this case radiation, is the external force that defines the matter property. Besides, most temperature models have a minimum time scale of daily estimations and the conceptual limitation that energy available for melt is not linearly related to positive air temperatures (Hock and Holmgren, 1996). Therefore, they are not able to reproduce daily pattern fluctuations which are important in tropical regions. Moreover, Kuhn (1987) performed a theoretical analysis of energy budget and melting conditions, and showed that melting may happen at air temperatures as low as −10 • C or as high as 10 • C. Hence, radiation must be included in melting models.
Energy balance melt models are based upon the assumption that at freezing temperature any surplus of energy at the surface air interface will be used for melting, and the energy available for melting is then related to the latent heat of fusion. The energy balance method has a strong physical background and was successfully applied in different studies (Molg et al., 2008;Cullen et al., 2007;Molg and Hardy, 2004;Wagnon et al., 1999). Some studies use results from energy balance models to calibrate temperature models (Carenzo et al., 2009). Basic energy balance models were applied to simulate snowmelt in Nordic glaciers  in the Alps (Sicart et al., 2008), New Zealand (Anderson et al., 2010) and in the Andes (Sicart et al., 2011). There are also more advanced energy models that divide the snowpack into layers and then apply a 1-D mass and energy balance to each layer in order to predict temperature profiles. The 1-D energy balance model SNTHERM was applied to analyse the variation of soil temperature with snow cover and improving roads maintenance (Fu et al., 2009). The land surface model ISBA-ES was coupled to models SAFRAN (meteorological model) and MODCOU (hydrogeological model) to simulate spring and summer flows in the French Alps (Lafaysse et al., 2011). However, the above mentioned models face the main limitation of detailed data requirements which are difficult to obtain, and sometimes it can be obtained only for limited periods of time. Long-wave radiation (LWR) is the energy term most difficult to obtain. Although it can be measured using a pyrgeometer, this is a very expensive instrument rarely available at meteorological stations (Dos Santos et al., 2011). LWR provides large amounts of melt energy for high-albedo snow surfaces and may even dominate in the energy balance for under cloudy skies (Sicart et al., 2010). Such fact has a special importance for the study area, as the local summer coincides with the wet season characterised by cloudy skyes and LWR controls the melting seasonal changes. Although simple equations were developed for its estimation considering vapour pressure and temperature (Kruk et al., 2010), those equations have empirical constants calibrated for specific locations. Hence, new alternatives are needed to estimate either energy fluxes or the whole energy balance using more accessible data.
ANN are mathematical structures able to represent complex non-linear relationships between input and output by imitating functioning of neurons in a human brain. In the last years, ANN have been successfully applied in hydrological studies. They were used for sediment studies (Kisi et al., 2012), weather forecast downscaling (Hoai et al., 2011), rainfall forecasting (Hung et al., 2009), river flow estimations (Akhtar et al., 2009;Shamseldin, 2010;Huo et al., 2012), litoral drift predictions (Singh et al., 2008) and evapotranspiration (Cobaner, 2011;Dai et al., 2009). Although widely applied to hydrology, almost no studies applied ANN to snow and glaciated areas. Yilmaz et al. (2011) applied ANN to estimate flow in a snow dominated mountainous basin in Turkey, but its time step was limited to daily scale and the model was not able to reproduce the yearly pattern, hence the model was developed as a seasonal model.
The present study developed different ANN models able to represent the nonlinear relations between common meteorological parameters, e.g. temperature, short-wave radiation or relative humidity, and other energy fluxes and the energy balance for a given time. Such models allow estimating energy balance and the energy available for melt. Relating the energy available for melt with the latent heat of fusion easily provides the potential SGM rate per unit of area. First, ANN models were developed for the estimation of energy fluxes difficult to obtain and usually unavailable, e.g. LWR sensible heat flux and latent heat flux. Then, ANN models were developed for estimating the whole EB within one single step. The ANN models were trained using measured data from the Bolivian glacier Zongo located in the outer tropics. Then, they were validated with measured data from the Equatorian glacier Antizana located in the inner tropics. The present research also developed MLR models for the estimation of energy fluxes. The performance of ANN models for estimating energy fluxes was compared against the MLR models. Both the ANN and MLR models were validated with data measured from Alpine glaciers. Finally, the validated models were applied to the Bolivian glacier Condoriri. The SGM estimation was compared with estimations from the MLR model and against enhanced temperature index models. This research is not only among the first ones to estimate SGM at short time step without complex data, but also among the first ones to implement ANN technologies in tropical glaciers in a developing country. The results will allow for easily predicting future SGM at any time. One main contribution of the present study is that it will allow overcoming the problem of data scarcity by proposing several models to estimate energy balance under different input data scenarios.

Study area
The present study considered three tropical glaciers: Antizana, Zongo and Condoriri. The Condoriri and Zongo glaciers are both in the outer tropics in the Bolivian Andes at some 13 km from each other, while the Equatorian glacier Antizana is located in the inner tropics at some 2100 km from Zongo and Condoriri. In order to test the validity of the model to glaciers out the tropical range, additional data from the Alpine glaciers Argentière and Saint-Sorlin were also used ( Fig. 1).
It is important to stress that the inner tropics and the outer tropics have different climatological behaviour. The humid inner tropics are characterised by a thermically homogenous atmosphere with slightly variable humidity throughout the year; temperature and humidity variations are not large enough to characterise a pronounced seasonal regime (Favier et al., 2004). On the other hand, the outer tropics have two marked seasons: one wet season characterised by precipitation and cloud cover coincident with the austral summer (November-March), and one dry season with clear skies during winter (April-September) (Sicart et al., 2005). The shortwave radiation (SWR) in Zongo and Condoriri has maximum values in winter due to the clear skies, while during summer it has lower values due to the cloudiness that attenuate it. On the other side, long-wave radiation has lower values in winter and higher values in summer due to the radiation reflected by clouds. The temperature reaches its highest values in summer and lowest values in winter; but due to the high altitude it is common to have frozen temperatures even in summer.
The Zongo glacier with altitudes ranging from 4900 to 6000 m a.s.l. (above sea level) is located in the Huayna massif (16 • 16 S, 68 • 10 W) and is part of a 3.7 km 2 basin with the main limnimetric station at 4830 m a.s.l. (Sicart et al., 2007). It is located some 30 km north of La Paz. The Antizana glacier with altitudes between 4800 to 5700 m a.s.l. is located 40 km east of Quito, on the NW slope of the Antizana volcano (Favier et al., 2004). Both glaciers are being monitored since 2003 within the project GLACIOCLIM (http://www-lgge.ujf-grenoble.fr/ServiceObs/).
The Condoriri glacier with altitudes from 4400 to 5200 m a.s.l. (16 • 11 S, 68 • 13 W) has the shape of a condor with open wings and provides water for the cities of El Alto and La Paz. Condoriri glacier along with Huayna and Tuni glaciers are currently being studied under the GRANDE project that will allow researchers to comprehend what is happening to the glaciers and to predict future scenarios (http://grande.civil.tohoku.ac.jp/indexe.html). One major problem of the above mentioned glaciers is the lack of data. In July 2011 GRANDE project installed weather stations around the mentioned glaciers. Although the stations will provide current data at small time intervals, the measured parameters are not enough to perform a complete energy budget.
Saint-Sorlin (45 • 10 N, 6 • 10 E) and Argentière (45 • 95 N, 6 • 98 E) glaciers are located in the Western Alps of France. They are monitored also within the project GLACIOCLIM. Meteorological stations have been set up close to both glaciers in order to study the relationships between climate change and mass balance fluctuations (Six et al., 2009). While Antizana Zongo and Condoriri are located in tropical latitudes in the Southern Hemisphere, Saint-Sorlin and Argentière are located at septentrional latitudes, with lower elevation, marked seasonality and well-defined periods of ablation and accumulation.

Data
Before applying the models to the Condoriri glacier, they were validated against data measured at the Equatorian glacier Antizana located in the inner tropics. Data for the glaciers Zongo and Antizana was obtained from meteorological stations installed and monitored by the French project GLACIOCLIM. Those stations have a data logger Campbell-ORE23x. They record every 30 min several meteorological parameters, such as short-wave radiation, long-wave radiation, temperature, relative humidity, wind speed and wind direction. The complete database consists of nearly 78 000 data points covering the years 2003-2009. The consideration of the number of data used for training the ANN models is an important factor. Too large ANN tend to overfit, while ANN with too few data do not contain enough processing elements to correctly model the input data set and tend to underfitting the data; both of these situations result in poor generalisation (Barnard and Wessels, 1992). Longer data sets should provide better ANN models, since it provides additional knowledge, thus considering more possibilities when training the model (Hertz et al., 1991). Hence, it is important to have at least a minimum number of training data. The present study considered that the number of training data should be at least 10 times the number of ANN weights (Baum and Hausler, 1989). It is important to consider that extracting more information from data does not always improve the models, indicating the need for improved data and models for critical times. Neal et al. (2011) suggest developing different models for day and night. The present study focused on daytime hours, defined as the ones when incoming short-wave radiation (ISWR) is higher than 20 W m −2 (Hu et al., 2012). Another limitation of ANN models is that the limits of the training data may impose a limitation to its application towards data outside the training data limits (Hettiarachi et al., 2005). In the present study the temperature (t) and relative humidity (RH) of the training data have lower limits than the validation data (Table 1). The validation data has an average RH 20 % higher than the training data. Although both training and validation data have similar average temperature, the instant maximum of the validation data is 5 • C higher. Nevertheless, the training data contains high values of the main source of energy, i.e. solar radiation, as high as the assumed solar constant of 1366 W m −2 (Spokas and Forcella, 2006). Data for the Condoriri glacier was obtained from the Condoriri weather station installed by the GRANDE project. This station installed in July 2011 has a data logger HOBO-U30 that records every 10 min several meteorological parameters like short-wave solar radiation, wind velocity, relative humidity, temperature and rain. The present study used more than 61 300 data points from July 2011 to August 2012, covering one hydrologic year.
Data for the additional validation glaciers of Argentière and Saint-Sorlin was obtained from meteorological stations installed and monitored by the French project GLACIO-CLIM. These stations provide SWR, LWR, relative humidity, temperature, wind speed and net radiation.

Energy model
The estimation of SGM by energy methods is based on the assumption that once the glacier reached freezing temperature any surplus of energy is used for melting . The method may be applied either at single where Q M is the energy flux available for melting, SW in is the incoming short-wave radiation flux, LW in is the incoming long-wave radiation flux, SW out is the outgoing shortwave radiation flux, LW out is the outgoing long-wave radiation flux, Q H is the sensible heat flux, Q L is the latent heat flux, and Q O are other minor heat fluxes like rain. Usually short-wave radiation is measured at the site. Sometimes it may be estimated by a valid relation considering the location (latitude and longitude), the Julian day, possible cloudiness and time exposed to solar radiation which is influenced by the local topography. In the present study, both LWR and SWR were obtained from Zongo station.
Sensible heat is calculated as function of the wind speed and temperature ) (Eq. 2): where C p is the specific heat air at constant pressure, k is the Von Karman constant, P is the atmospheric pressure, u is the wind speed, P o is the standard atmospheric pressure, Z is the instrument height, Z ow is the roughness for wind logarithmic profile, Z ot is the roughness for temperature logarithmic profile, T f is the freezing temperature, P is the air density. Latent heat is calculated as function of the wind speed and humidity where L is the latent heat flux of evaporation, e 2 is the vapour pressure at 2 m, e o is the vapour pressure at melting surface, Z oe is the roughness parameter for vapour pressure logarithmic profile. Vapour pressure was estimated as function of relative humidity (Allen et al., 1998) where e S is the saturation vapour pressure and RH is the relative humidity [%] (Eq. 5).
where T is temperature in Celsius degrees. Stability correction described by the Richardson number (Eq. 6) was applied to the turbulent fluxes (Sicart et al., 2005;Favier et al., 2004).
where T and u are absolute air temperature [K] and wind speed [m s −1 ] at the measurement level Z, g is the acceleration of gravity (9.79 m s −2 ), T s is the surface temperature and Z om is the roughness length [m]. The non-dimensional stability functions for momentum ( m ), heat ( h ), and moisture ( v ) were expressed in terms of Ri. Ri positive (stable) Ri negative (unstable) The sensible heat of fusion was multiplied by ( m h ) −1 and the latent heat of fusion was multiplied by ( m v ) −1 . The other energy fluxes were neglected, since they represent a very small percentage of the total. Then, the energy available for melt was converted into its water equivalent by relating to the water latent heat of fusion (Eq. 9).
where WE is the water equivalent per square meter [mm s −1 m −2 ], L f is the water latent heat of fusion assumed as 334 J g −1 . Hock and Holmgren (1996) suggested roughness values for glacier areas of Z ow = 0.0027 m, Z ot = 0.000027 m, Z oe = 0.000027 m. Sicart et al. (2011) studied the Zongo glacier and suggest a value of Z ow ranging from 1 to 10 mm and Z ot = Z oe = Z ow /100, while other studies consider them as calibrating variables. Sicart et al. (2005) found that both latent and sensible heat fluxes in Zongo are small and they play a minor role in the total energy balance since radiation supplies most of the melting energy. Also Van As (2011) found that solar radiation is the main source of melting energy, and the errors of assuming constant roughness are negligible. Thus, it can be assumed that the uncertainties of considering constant roughness values are small and without much influence.
The energy balance for each time step was estimated using data from the glaciers Zongo and Antizana. Then, the EAM was estimated by comparing with previous time step energy in order to consider refreezing effects; in case the energy balance was negative, it was assumed as freezing  and the next time step energy must compensate such freezing before allowing for melting. Once there was enough energy for melting, such energy was related to the latent heat of fusion for water which was assumed 334 J g −1 in order to get the melting water equivalent for that time step. Ideally, snowmelt and ice melt should be differentiated, since they have slightly different L f depending on the amount of liquid water in the snow (Sing and Vijay, 2001;Guttman, 1907). However, since the differences are too low, it was assumed that both snow and ice have the same L f .

Artificial neural networks
ANN are approximation methods that imitate the functioning of the human's brain. The brain may be idealised as a highly complex non-linear and parallel computer with the capability to perform computations by organizing its neurons and building up its own rules through learning process. In analogy, ANN may reproduce multi-variable functions by arranging processing elements (neurons) interconnected according to certain rules that may change in order to find the optimal ones (learning). The most popular type of ANN is the multilayer perceptron (MLP), which is a feedforward network that has interconnected nodes (neurons) arranged into three layers: input layer, a hidden layer and an output layer (Fig. 2) (de Vos, 2013;Hung et al., 2009). The input layer sends the input vector X of signals x i to the hidden layer. The hidden layer enables the network to learn by extracting meaningful features from the input. Each neuron processes its output y j by summing its input signal x i multiplied by its respective weight w ij and a given threshold a o (Eq. 10).
The output of each neuron may go to the next hidden layer or to the output node (if there is only one hidden layer). The main characteristics of a MLP are (a) that each neuron includes a soft nonlinearity (sigmoidal logistic function which is described in Eq. (11), (b) its layered architecture allows to learn by progressively extracting information from the input, and (c) its high degree of connectivity so that one element of a given layer feeds all the nodes of the next layer.
In this study the MLP was trained with the back-propagation algorithm of the software Waikato Environment for Knowledge Analysis (WEKA) version 3.6.6 (Hall et al., 2009). The performance was evaluated by the non-overlapping test set selection cross-validation method, also known as the k fold method. This is one of the most popular validation methods, and can be described in the following five steps (Gascon-Moreno et al., 2012): 1. The total available data n is divided into k nonoverlapping data subsets C1, C2, ..., Ck, also known as folds.
2. One fold is used for validation, while the remaining folds are used as training data.
3. The created model is tested with the testing fold. This test generates an 4. Steps two and three are repeated k times, so that every fold is used once as validation set.
Hydrol. Earth Syst. Sci., 17, 1265-1280, 2013 www.hydrol-earth-syst-sci.net/17/1265/2013/ While normal hold out testing methods may produce biased results due to the data partitioning, this method gives a much fair and unbiased estimation (Bengio and Grandvelet, 2004). The popularity of the method grew so fast, that studies also compared its performance by using different number of folds. While early ideas suggested using 10 folds, it was found that the differences between 5 and 10 folds are not significant (Uguz and Kodaz, 2011;Iliadis et al., 2011). Markatou et al. (2005) stated that 4 (folds) is a reasonable number that provides fair estimations, hence the present study used 4 folds. The performance of each model was evaluated in terms of correlation coefficient (correlation), mean absolute error MAE (Eq. 13) and root mean squared error (RMSE) (Eq. 14).
where P i is the predicted value, T i is the target value, n is the total number of samples.
Since SWR may easily be measured by using a pyrometer, which is a simple instrument usually available in meteorological station, e.g. the main study basin of Condoriri is equipped with pyrometers, SWR was considered as an input data. The energy flux not available is long-wave radiation, which is rarely available due to the high cost and complexity of pyrgeometers (Choi et al., 2008). Hence, as a first step, ANN models were developed to simulate incoming LWR (ILWR) and global LWR (GLWR). Since all the ILWR models consider temperature and vapour pressure (Dos Santos et al., 2011), such data was considered as input data. Temperature is measured, and vapour pressure is estimated from relative humidity, as previously explained. Different models were developed considering temperature in Celsius degrees, in Kelvin degrees and to the 4th power. Besides, it was important to consider that ILWR has to be corrected considering the cloud effect. Since previous studies showed a correlation between relative humidity and cloud cover (Walcek, 1994), it may be assumed that including relative humidity may represent the cloudiness effect. In total 6 ILWR models and 6 LWR models were developed, as showed in Table 2 groups a and b.
Since latent heat flux and sensible heat flux are estimated using temperature (t), relative humidity (RH) and wind speed (u) , in order to perform additional test of the methodology it was decided to develop additional ANN models for the estimation of LHF and SHF considering t, RH and u, as showed in Table 2 group c.
Although the independent models allows to estimate the energy balance in two steps, i.e. estimate each energy flux and then estimate the energy balance, additional ANN models were developed in order to estimate the whole energy Tc is the temperature in Celsius degrees. Tk is the absolute temperature in Kelvin degrees. RH is the relative humidity. u is the wind speed. The data with the subscript "1" denotes data from a previous time step. The superscript 4 means to the 4th power. ILWR is incoming long-wave radiation. GLWR is global long-wave radiation. SHF is sensible heat flux. LHF is the latent heat flux. EB is energy balance. EAV is energy available for melt.
budget and the EAM within one single step. Table 2 shows the developed models and data required for each one. The predictive performance of the ANN models were compared against predictive models based on the MLR technique, which is a recognised statistical technique widely used in engineering and science problems dealing with exploring the relationship between two or more variables (Agha and Alnahhal, 2012;Riad et al., 2004). The MLR were developed using the analysis tool pack of MS Excel 2010.
Finally, the models with the best performance were used for estimating SGM. SGM at Condoriri not only was estimated with both the ANN and MLR models, but also with an enhanced temperature index model that considers temperature and short-wave radiation (Eq. 15) (Pellicioti et al., 2005).
If T > T t : The models with the subscript "a" were developed also considering data from a previous time step. The term techniques refers to whether the model used was ANN or MLR. The term criteria refers to the criteria used for evaluation of the performance.
Otherwise M is 0. Where M is the melting and α is albedo. TF and SRF are empirical coefficients and T t is the threshold temperature for melting. TF and SRF were assumed 0.05 mm h −1 • C −1 and 0.0094 mm W −1 h −1 as suggested by Pelliccioti et al. (2005). Since melting may happen at different temperatures within a range of −10 to 10 • C (Kuhn, 1987), a sensitivity analysis was performed considering different melting temperatures within the mentioned range. When estimating SGM via the energy models the conversion from EAM into water equivalent was allowed if the temperature of the given time step was higher than the threshold. It is important to note that the estimated SGM is the water equivalent. Thus, it was converted into mm of ice by assuming a mean surface snow glacier density of 0.4 g cm −3 (Ginot, 2001)

Results and discussion
Tables 3 and 4 show the comparison of performance between ANN and MLR for glaciers Zongo and Antizana, respectively. ANN and MLR have similar performance when estimating LWR at Zongo. For the incoming LWR, the MLR The models with the subscript "a" were developed also considering data from a previous time step. The term technique refers to whether the model used was ANN or MLR. The term criteria refers to the criteria used for evaluation of the performance. The term cor refers to correlation. models have slightly better performance. The ANN models have better performance for the global LWR. When validating the incoming LWR models at Antizana, the MLR models have slightly lower errors, but they have almost the same correlation; in the first three models the MLR correlation is 0.01 higher, and in the other 3 models the ANN correlation is 0.01 higher. It was found that there is no difference between using relative or absolute temperature, i.e. temperature may be used as Celsius or Kelvin degrees. The use of temperature to the fourth power does not provide significant differences. This is an advantage compared with empirical equations that require the conversion of temperature to kelvin degrees and its elevation to the fourth power. When applying the model to the Condoriri location it can be seen that it clearly reflects the seasonality (Fig. 3). The winter period between May and August has the lowest LWR fluxes. During the spring period between September and November the LWR fluxes increase, and during the wet season between December and March the LWR fluxes reach the maximum values. Analysing the models developed for representing the sensible heat flux, it can be noted that the ANN has a much better performance than the MLR in both Zongo and Antizana. The ANN reaches correlations of 0.98 and 0.85 for Zongo Fig. 3. Monthly incoming long-wave radiation flux estimated at Condoriri. The white rectangle represents the inter quantile range (IQR). The thick line inside the IQR represents the mean, and the circles represent the outliers defined as extreme events with a probability lower than 1 %. and Antizana, respectively, while the MLR has correlations of 0.86 and 0.73 for Zongo and Antizana, respectively. Both the ANN and the MLR have similar magnitude errors. The case of the latent heat flux is different. Although the MLR has a better correlation, it also has higher errors both in Zongo and Antizana. It may be assumed that the MLR tends to over fit its training data. Hence, it may be assumed that the ANN provides a better generalisation of both turbulent fluxes.
Analysing the models for representing the whole EB, SWR is the most important term and neglecting it gives correlations lower than 0.56. Combining only short-wave radiation and temperature leads to coefficients of 0.85 and 0.92 for the ISWR and GSWR, respectively. Including relative humidity improves the correlation by about 10 %, reaching correlations of 0.92 and 0.98 for ILWR and GLWR, respectively. Including relative humidity also decreases the errors. The inclusion of wind speed does not provide relevant improvement; moreover, in the validation data set of training the inclusion of wind speed slightly reduces the correlation. Both the ANN and the MLR have the same correlations and similar errors. As in the previous models, the ANN models have lower errors than the MLR models when applied to the validation data set.
Analysing the energy available for melt, i.e. considering the heat of previous freezing (Models EAV7 to EAV10 a ), the training models considering only current data (EAV7 to EAV10) have correlations lower than 0.88. Models that also consider data from a previous time step (EAV7 a to EAV10 a ) have higher correlation and lower errors. The models that neglect SWR have correlation lower than 0.51. Hence, energy available for melt may not be accurately estimated without SWR. As in previous models, the ANN models perform better in the validation data set with lower errors, thus providing a better generalisation of the phenomena. Considering models with previous time steps (models EAV7 a to EAV10 a ), ANN models also perform slightly better than MLR, e.g. ANN model reaches a correlation 0.9 while MLR reaches correlation of 0.88. In the models neglecting SWR, MLR performs better than ANN both in training and validation data; nevertheless, it still has correlation lower than 0.63. Applying the ANN models to a different location with data outside the training set limits provides results with similar correlation, which may be considered as a hint that the ANN provides a good estimation of the energy fluxes and energy balance; however, they provide higher errors and an underestimation of total energy. MLR also underestimated the total energy with higher errors. All these results indicate that ANN provides a better generalisation and that they are less sensitive to errors when applied to data outside its training limits. It is important to note that the magnitude of the errors and the underestimation is related to the relevance of the data outside of the training range. For instance, the models for LWR are driven mainly by temperature and relative humidity, where Antizana input data is higher that the training data; therefore, the application of those models to Antizana provides results with the highest errors. On the other side, models where SWR is the most important input data have results with lower errors. Instead of developing two different models for each location or applying a correction coefficient, it was decided to improve the existing models by including the data from Antizana, hence incorporating its knowledge into the ANN. Another approach suggested and tested in the present study for dealing with extrapolation data is the prior transformation of data to its logarithmic (Hettiarachi et al., 2005). Table 5 shows a comparison of ANNs developed with the prior logarithmic transformation with ANN and MLR models developed using data from Zongo and Antizana. These techniques were applied for developing a new ILWR4 model. That model was selected as it had the lowest performance when validated at a different location. The ANN model that used data from both Zongo and Antizana reached a correlation of 0.82. The MLR model with data from Zongo and Antizana reached correlation of 0.81 and higher errors.

V. Moya Quiroga et al.: Snow glacier melt estimation in tropical Andean glaciers
Although the ANN using the logarithmic transformation has higher correlation and lower errors in the training set, it does not provide relevant improvements when applied to Antizana. It is important to consider that the improved ANN model that used data from Zongo and Antizana was developed with data from tropical latitudes and located at high elevation; thus, it considers a wide range of input data including extreme solar radiation of 1349 W m −2 which is near the solar constant (1360 W m −2 ) and temperatures up to 16 • C. Those values are much higher than normal values from high latitude glaciers like the Alps. In order to validate the new developed models, an additional comparison was performed against measured data from the French glaciers Argentière and Saint-Sorlin (Table 6). This validation reached correlation higher than the training data and similar magnitude errors. The new ANN models have a better performance than the new MLR models. Besides, the ANN models have the same performance with similar magnitude errors for both locations, while the MLR models have a better performance when applied to Argentière glacier. Therefore, it may be assumed that ANN models provide a better generalisation than MLR models.
SGM was then estimated at the tropical glaciers Zongo and Antizana using different ANN models, the energy balance and the enhanced temperature method. This comparison also considered different melting temperatures between −10 and 10 • C. Considering that energy balance was successfully applied in previous studies (Molg et al., 2004(Molg et al., , 2008Cullen et al., 2007;Wagnon et al., 1999) and its results are used to calibrate temperature based equations (Carenzo et al., 2009), it was decided to use results from the energy balance as reference SGM. It was found that in both locations the temperature model overestimated the total SGM (Tables 7 and 8).
On the other hand, the ANN models provide results closer to the ones based on the energy balance. All the models have the same trend with nearly constant SGM for melting temperatures lower than 0 • C and an abrupt decrease for higher temperatures. The ANN models have same performance and same magnitude errors in both locations; on the other hand, the enhanced temperature method has better performance at the Zongo glacier. Thus, it may be assumed that the ANN models provide a better generalisation and may be used at different locations, while the enhanced temperature method has to be checked and calibrated for each location.
Comparing the SGM estimated at Condoriri with the different ANN models that can be used with the Condoriri data and considering different melting temperatures (Table 9), it may be noticed that all the models have the same trend according to the melting temperatures. The models that neglect solar radiation provide the highest predictions, while the models that neglect temperature provide the lowest predictions. The other nine models provide predictions that have values within a confidence interval limited by the upper confidence limit (UCL) and the lower confidence limit (LCL) as described by Eqs. (16) and (17).
where SGM av is the average SGM, and SD is the standard deviation.
Hence, it may be assumed that solar radiation plays an important factor and its neglecting may lead to higher errors. The model EAV3 is the closest model to the average SGM.
Total melting was estimated at Condoriri with 3 SGM models: the ANN with highest correlation, the MLR with highest correlation and the enhanced temperature method of Pelliccioti. It may be noticed that the temperature method overestimates the melting about 33 %. Hence, this model was adjusted by including a correction coefficient of 0.75 (Fig. 4). The ANN and the MLR provide similar values with total differences of about 0.3 m. Such results were obtained considering a fixed melting threshold temperature of 0 • C. Considering different melting temperatures between −10 and 10 • C it may be noticed that all models have the same trend for different melting temperatures. Considering melting threshold temperatures between −10 and 5 • C, the total SGM is about 6.45 m with differences of about 0.4 m. Melting threshold temperatures higher than 5 • C highly decrease the total melting. Assuming a melting threshold temperature of 10 • C would produce a total melting of 1.80 m. Comparison of the SGM considering freezing effect and without considering it results in a total yearly difference of 0.27 mm WE. Figure 5 shows the hourly SGM rates for the months of March, June, September and December in statistical boxplots (Frigge et al., 1989). Those months were selected as representative for each season, since they are the ones when solstice and equinox occur. These plots show the seasonality. The most important change is not the highest melting rate, but the number of hours that the location is exposed to melting. During the winter month of June, the melting times are limited from 08:00 to 14:00 LT (Bolivian local time). During the summer month of December the melting hours are from 06:00 to 17:00 LT which represents five additional melting hours. The months of March (fall) and September (spring) have the same hours of exposure to melting (from 07:00 to 16:00 LT) and similar melting rates.

Conclusions
The present research developed artificial neural networks (ANN) models to estimate energy fluxes, the energy balance and the energy available for melting by the use of simple and easy to obtain data (short-wave radiation, temperature and relative humidity). Different models were developed considering different input data. Although ANN models and multiple linear regression (MLR) models have similar performance, the ANN provide slightly better performance when applied to different geographical locations with data outside the training set range. An additional benefit of ANN and MLR models compared with physical equations for estimating other energy fluxes, e.g. long-wave radiation, is that they do not require previous conversion of measured temperature into absolute temperature. Moreover, ANN and MLR do not require empirical coefficients to be calibrated for different locations.
Although the model with the highest number of variables is the most accurate, the other models also have a good performance with slightly lower performance. This presents a special advantage, since the proposed models allow estimating SGM under different data availability scenarios.
Short-wave radiation is the main source of energy and the data with the most influence on the predictive performance of the ANN and MLR models.
The limits used in the training set have a strong influence on the performance of the models. Application of ANN models to data outside the limits of the training set reduces the accuracy of the predictions. It is important to note that the magnitude of the errors and the underestimation is related to the relevance of the data outside of the limits range. MLR are more sensitive to the limitations of the training set and present higher errors when applied to different locations with data higher than the training set. , winter (June), spring (September) and summer (December). The white rectangle represents the inter quantile range (IQR). The thick line inside the IQR represents the mean, and the circles represent the outliers defined as extreme events with a probability lower than 1 %. It has been shown here that increasing the range of the training set, i.e. the addition of a wider domain of knowledge, provides more improvement in the predictive capacity of the ANN than other techniques like prior transformation of training data.
The training set used for developing the ANN includes extreme values of SWR and temperature. The SWR includes values as high as the solar constant. Besides, it used data from the humid inner tropics and the arid Andean plateau. It considers different relative humidity values and temperatures up to 16 • C. The ANN models were successfully tested at a different geographical location (the French Alps).
Temperature methods for estimating SGM require the calibration of coefficients. The use of uncalibrated coefficients may induce errors higher than 30 %. On the other hand, estimation of SGM via energy balance (either ANN or MLR) does not require further calibration.
The present paper presents a novel methodology for estimating the complex process of SGM at different hours in tropical glaciers using easy to obtain data. The developed models reflect the daily pattern of SGM which is very important for the study of tropical glaciers, since the daily variation is greater than the yearly one.