Skip to main content Accessibility help
×
Home
Hostname: page-component-7ccbd9845f-w45k2 Total loading time: 1.61 Render date: 2023-01-31T17:26:55.594Z Has data issue: true Feature Flags: { "useRatesEcommerce": false } hasContentIssue true

Mass balance and area changes of glaciers in the Cordillera Real and Tres Cruces, Bolivia, between 2000 and 2016

Published online by Cambridge University Press:  12 December 2019

Thorsten Seehaus*
Affiliation:
Institute of Geography, Friedrich-Alexander-University Erlangen-Nuremberg, Wetterkreuz 15, 91058Erlangen, Germany
Philipp Malz
Affiliation:
Institute of Geography, Friedrich-Alexander-University Erlangen-Nuremberg, Wetterkreuz 15, 91058Erlangen, Germany
Christian Sommer
Affiliation:
Institute of Geography, Friedrich-Alexander-University Erlangen-Nuremberg, Wetterkreuz 15, 91058Erlangen, Germany
Alvaro Soruco
Affiliation:
Instituto de Investigaciones Geologicas y del Medio Ambiente, Universidad Mayor de San Andrés, Calle 27 Cota Cota, La Paz, Bolivia
Antoine Rabatel
Affiliation:
CNRS, IRD, Grenoble INP, Institut des Géosciences de l'Environnement (IGE, UMR5001), Univ. Grenoble Alpes, 38000Grenoble, France
Matthias Braun
Affiliation:
Institute of Geography, Friedrich-Alexander-University Erlangen-Nuremberg, Wetterkreuz 15, 91058Erlangen, Germany
*
Author for correspondence: Thorsten Seehaus, E-mail: thorsten.seehaus@fau.de
Rights & Permissions[Opens in a new window]

Abstract

Climate change has led to a significant shrinkage of glaciers in the Tropical Andes during the last decades. Recent multi-temporal quantifications of ice mass loss at mountain range to regional scale are missing. However, this is fundamental information for future water resource planning and glacier change projections. In this study, we measure temporally consistent glacier area changes and geodetic mass balances throughout the Bolivian Cordillera Real and Tres Cruces based on multi-sensor remote-sensing data. By analyzing multi-spectral satellite images and interferometric SAR data, a glacier recession of 81 ± 18 km2 (29%; 5.1 ± 1.1 km2 a−1), a geodetic mass balance of −403 ± 78 kg m−2 a−1 and a total ice mass loss of 1.8 ± 0.5 Gt is derived for 2000–2016. In the period 2013–2016, ice mass loss was 21% above the average rate. A retreat rate of 15 ± 5 km2 a−1 and a mass budget of −487 ± 349 kg m−2 a−1 are found in this more recent period. These higher change rates can be attributed to the strong El Niño event in 2015/16. The analyses of individual glacier changes and topographic variables confirmed the dependency of the mass budget and glacier recession on glacier aspect and median elevation.

Type
Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright
Copyright © The Author(s) 2019

1. Introduction

Glaciers and ice caps outside of the polar ice sheets are defined as essential climate variables by the Global Climate Observing System (GCOS) and key indicators for climate change by the Intergovernmental Panel on Climate Change (IPCC, 2014). In particular, glaciers in low latitudes, also known as tropical glaciers, are very sensitive to changes in climatic conditions (e.g. Kaser and Osmaston, Reference Kaser and Osmaston2002; Rabatel and others, Reference Rabatel, Machaca, Francou and Vincent2006, Reference Rabatel2013). Small changes in temperature and precipitation strongly influence the mass balance of tropical glaciers since they are close to melting conditions throughout the year (Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003; Vuille and others, Reference Vuille2008).

The El Niño Southern Oscillation (ENSO) is known to affect the precipitation pattern in the Tropical Andes. The decadal precipitation variability appears larger than the long-term trend, which is additionally regionally contrasted. Indeed, in Ecuador and northern Peru, a long-term increase in precipitation since the mid-20th century is documented (Vuille and others, Reference Vuille, Bradley, Werner and Keimig2003). In contrast, a decreasing trend is reported for the Bolivian Altiplano region (Vuille and others, Reference Vuille, Bradley, Werner and Keimig2003; Haylock and others, Reference Haylock2006). Moreover, between 1939 and 2006, Vuille and others (Reference Vuille2008) documented a near-surface air temperature increase of 0.68 °C for the Tropical Andes (1°N to 23°S), which agrees with local observations by other studies (e.g. Vuille and Bradley, Reference Vuille and Bradley2000; Mark and Seltzer, Reference Mark and Seltzer2005). Gilbert and others (Reference Gilbert, Wagnon, Vincent, Ginot and Funk2010) analyzed englacial temperature measurements in a 138 m deep borehole close to the summit of Illimani (6340 m a.s.l., Bolivia). The authors derived a general warming trend of 1.1 ± 0.2 °C for the 20th century, with a higher rate in recent years (4.3 ± 1.4 °C per century for 1985–1999). As a consequence of the changes in climatic conditions, the equilibrium line altitude (ELA) of the glaciers in this region rose (Rabatel and others, Reference Rabatel2013; Réveillet and others, Reference Réveillet, Rabatel, Gillet-Chaulet and Soruco2015). At the Quelccaya Ice Cap, Peru, an ELA shift from 5275 to 5414 m a.s.l. between 1975 and 2012 was observed by Hanshaw and Bookhagen (Reference Hanshaw and Bookhagen2014), and at Zongo Glacier, Bolivia, a temperature sensitivity of the ELA of 150 ± 30 m °C−1 was estimated (Lejeune, Reference Lejeune2009).

In the Tropical Andes, glaciers and ice caps are an important water resource and regulate the seasonal runoff, particularly in central and southern Peru and Bolivia (Vuille and others, Reference Vuille2018). No seasonal snow cover exists outside the glaciers, which could act as a buffer for the dry season runoff (Lejeune and others, Reference Lejeune2007). During the dry season, from April until September, the glacier meltwater significantly contributes to the drinking water supply of Bolivia's capital La Paz (e.g. Soruco and others, Reference Soruco2015). Additionally, glacier runoff is highly important for irrigation, hydropower generation and water consumption of local communities (Vuille and others, Reference Vuille2018). The extreme drought in Bolivia in 2016 forced the administration to declare an emergency and to ration water (Schoolmeester and others, Reference Schoolmeester2018). The authors also reported a very high contribution of 61% of the glacier meltwater to the available water supply in La Paz during this extreme event. Buytaert and others (Reference Buytaert2017) estimated an even higher monthly maximum of 85.7% during drought years. On a multi-decadal timescale (1963–2006), Soruco and others (Reference Soruco2015) reported an average glacier meltwater contribution to the water supply of La Paz city, the administrative capital of Bolivia, and the neighboring city El Alto (together ~2.3 million inhabitants) of 27% during the dry season. The glacier recession in the Tropical Andes leads to temporarily increased meltwater runoff. However, on long-term scales, the runoff will decrease as glaciers retreat further (Vuille and others, Reference Vuille2018). Estimates by Huss and Hock (Reference Huss and Hock2018) indicate that the Bolivian glaciers have already reached or are close to ‘peak water’, the point in time with maximum meltwater runoff. Soruco and others (Reference Soruco2015) projected a reduction in runoff of 24% during the dry season if the glaciers disappear in the La Paz drainage basins.

The glacier area in the Bolivian Cordillera Oriental (Cordillera Apolobamba, Cordillera Real and Tres Cruces) reduced by 43% between 1986 and 2014 (Cook and others, Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) and even a higher retreat of ~50% was found in the small mountain range of Tres Cruces between 1975 and 2009 (Albert and others, Reference Albert, Kargel, Leonard, Bishop, Kääb and Raup2014). For the southern wet outer tropics, where the Bolivian Cordillera Oriental is located, a mass change rate of −0.44 ± 0.11 Gt a−1 (850 kg m−3 density scenario) was estimated by using radar remote-sensing data for the period 2000–2013 (Braun and others, Reference Braun2019). This estimate also includes largely glacierized areas in Peru and the authors do not provide information on detailed spatial scales since the study focused on a continent-wide mass-balance analysis. Detailed studies of glacier changes in Bolivia have been carried out at the Zongo, Chacaltaya and Charquini glaciers (e.g. Rabatel and others, Reference Rabatel, Machaca, Francou and Vincent2006; Soruco and others, Reference Soruco2009b; Réveillet and others, Reference Réveillet, Rabatel, Gillet-Chaulet and Soruco2015; Soruco and others, Reference Soruco2015). Volume changes between 1963 and 2006 of 21 glaciers in the Cordillera Real were measured by Soruco and others (Reference Soruco, Vincent, Francou and Gonzalez2009a) based on photogrammetric measurements. The authors estimated for nearly one-third of the Cordillera Real's glaciers (376 glaciers) an area loss of 48% and a volume loss of 43% throughout the study period using an empirical relationship between mass balance, exposure and average altitude based on the 21 monitored glaciers.

At Cordillera Real and Tres Cruces, no region-wide, spatially detailed and multi-temporal information on the contemporaneous glacier area and mass changes is available. Therefore, the scope of this paper is to monitor glacier evolution throughout the Bolivian Cordillera Real and Tres Cruces for the period 2000–2016 in order to obtain enhanced information on recent glacier changes. Our primary objectives are to:

  • Generate detailed glacier inventories for the years 2000, 2013 and 2016, which are temporally consistent to the mass-balance analyses

  • Measure regional and glacier-wide geodetic mass balances for the periods 2000–2013, 2000–2016 and 2013–2016

  • Analyze area change and mass-balance results in comparison with topographic and climatic variables

The time intervals are defined by the availability of the used remote-sensing data and complement the observation periods of other studies.

2. Study site

The Cordillera Real and the smaller mountain range Tres Cruces are home of the majority of the Bolivian glaciers. In agreement with Jordan (Reference Jordan1991), the glacierized regions of Mururata and Illimani were included in the Cordillera Real (Fig. 1). The analyzed mountain ranges extend in the NW-SE direction and have a length of ~200 km. They separate the wet Amazon Basin from the arid Altiplano Plateau. The study site is located in the outer tropics region as defined by Troll (Reference Troll1941), particularly in the southern wet outer tropics according to Sagredo and Lowell (Reference Sagredo and Lowell2012). Glacier surface albedo and the overall radiation budget, mainly controlled by cloud cover and solid precipitation, strongly force the surface mass balance of these glaciers (Wagnon and others, Reference Wagnon, Ribstein, Francou and Pouyaud1999; Sicart and others, Reference Sicart, Wagnon and Ribstein2005). During the wet summer season (January–April), accumulation, as well as snow and ice melt rates, are the highest, whereas in the dry winter season (May–August), less melting occurs and sublimation becomes important. The transition period (September–December) is characterized by sporadic precipitation events and high melt rates, since the solar radiation increases (Sicart and others, Reference Sicart, Hock, Ribstein, Litt and Ramirez2011; Rabatel and others, Reference Rabatel2012). ENSO events have also been identified as an important climate driver of the interannual variability of the glacier surface mass balance (e.g. Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003; Vuille and others, Reference Vuille2008; Rabatel and others, Reference Rabatel2013). El Niño events amplify the mass loss, whereas during La Niña events, lower mass loss or even positive mass balances are found (Wagnon and others, Reference Wagnon, Ribstein, Francou and Sicart2001; Vuille and others, Reference Vuille2008). Based on the Randolph Glacier Inventory (RGI) 6.0 (RGI Consortium, 2017), the glacierized area in the study region ranges between ~3700 (too low, caused by an artifact in RGI 6.0, see Section 5.1) and 6400 m a.s.l., with an extension of ~329 km2 (measured in UTM 18S, since the majority of the glaciers in RGI Region 16 are located in UTM 18S zone), corresponding to ~63% (according to RGI 6.0) of the glacierized area in Bolivia, consisting of ice caps, and mountain and valley glaciers. At three glaciers in the study region, continuous surface mass and energy-balance monitoring are (were) carried out (Zongo Glacier since 1991, Chacaltaya Glacier 1991–2008 and Charquini Sur Glacier since 2002) within the framework of a Franco-Bolivian collaboration through the French Service National d'Observation GLACIOCLIM. Zongo Glacier is a benchmark glacier of the World Glacier Monitoring Service (WGMS) and a reference Cryonet Station of the Global Cryosphere Watch led by the World Meteorological Organization.

Fig. 1. (a) Map of Bolivia. Dashed square indicates the study area, blue polygons show the glacier extent based on the RGI 6.0. Background: Natural Earth (b) Topographic map of analyzed mountain ranges and TanDEM-X (TDX) and Pléiades data coverage. Background: SRTM © NASA.

3. Data and methods

3.1 Area changes

A detailed comparison of the RGI outlines with high-resolution satellite imagery of the respective years (from Google Earth) revealed a misclassification of certain areas. As stated in the Technical Report of the RGI 6.0 (RGI Consortium, 2017), a more rigorous examination and delineation of the glacier outlines might reduce the mapped glacier area of RGI Region 16 by up to 20%. This discrepancy can bias area change computation when comparing the RGI outlines with other outlines.

Consequently, glacier outlines are mapped, using cloud-free LANDSAT images (Table S1, provided by the United States Geological Survey, USGS), for 2000, 2013 and 2016. To identify glacierized regions, the NDSI of top of atmosphere reflectance values and the band ratio (BR) (GLIMS Algorithm Working Group) of red and shortwave infrared bands of the digital number values are calculated. A combination of NDSI and BR (for deep shadows) information is used to identify the glacier outlines. By means of high-resolution satellite imagery (Google Earth) from the respective years, the threshold values for the NDSI of 0.8 and the BR of 1.7 for LANDSAT 5 TM data and 1.5 for LANDSAT 8 OLI data are carefully selected. The resulting binary masks are transferred into polygons and cropped into individual glacier basins using the ice flowshed delineations of the RGI. Afterward, the glacier inventories are then visually inspected and misclassified areas, such as proglacial lakes and snow patches, are removed by means of high-resolution satellite data (Google Earth). The surface area (S) of the final inventory and the individual glaciers is calculated in UTM 19S projection. Additionally, the topographic parameters (minimum, maximum and median elevation, slope and aspect) of the inventoried glaciers are calculated based on the void fill SRTM Version 3 digital elevation model (DEM) (NASA JPL, 2013, see also next section) according to RGI6.0 Technical Report (RGI Consortium, 2017).

3.2 Elevation changes and mass balance

In this study, multiple DEMs derived from Interferometric Synthetic Aperture Radar (InSAR) data are used to determine glacier surface elevation changes throughout 2000–2016.

The Shuttle Radar Topography Mission (SRTM) of the National Aeronautics and Space Administration (NASA) and the German Aerospace Center (DLR) acquired bi-static C-band InSAR data of landmasses between 60°N and 54°S from 11 to 22 February 2000. The goal of SRTM was to facilitate the generation of a nearly global, high-resolution and consistent DEM. Several DEM products were derived from SRTM data. In this study, the void-filled LP DAAC NASA Version 3 SRTM DEM (1 arcsec ground resolution) (NASA JPL, 2013), reprojected to UTM projection, is used. Further DEMs are generated using data sets from the TanDEM-X (TDX) mission of the DLR. The TDX mission has been acquiring bi-static X-band InSAR data since 2010 (Zink and others, Reference Zink, Bartusch and Miller2011). A nearly complete coverage of the study site by TDX data is available for 2013 and 2016. In order to reduce the biases on the elevation change calculation caused by seasonal surface elevation changes and SAR signal penetration (see Section 5.3), TDX data within the same season as the SRTM are selected for 2013. Unfortunately, in 2016, TDX data are only available for the period September to November. The potential bias caused by this seasonal offset is discussed in Section 5.3. The footprints of the TDX acquisitions used are shown in Figure 1 and the main characteristics are given in Table S1.

The TDX DEMs are processed according to Seehaus and others (Reference Seehaus2016) and Malz and others (Reference Malz2018). Along track acquisitions from the same path and date are concatenated and a differential interferogram is generated using the SRTM DEM as a reference DEM. Subsequently, the differential interferogram is filtered and unwrapped using either the branch cut or minimum cost flow algorithm (the best results are manually selected) and the differential phase is converted into differential height values. Finally, the reference DEM heights are added and the resulting DEM is geocoded. In the next step, the TDX DEMs are bi-linearly vertically co-registered to the SRTM DEM on stable ground. Stable points are defined as areas with a slope of <25° (using the SRTM DEM) and by masking out vegetation (NDVI threshold of 0.3) and water (NDWI threshold of 0.1) bodies by means of LANDSAT data. Afterwards, the TDX DEMs are horizontally co-registered to the SRTM DEM using the iterative algorithm of Nuth and Kääb (Reference Nuth and Kääb2011). To remove the remaining vertical offsets, the TDX DEMs are again bi-linearly vertically adjusted to the SRTM elevations and finally mosaicked to one single DEM, including a time stamp for the individual pixels.

The SRTM DEM and the resulting TDX DEM mosaics are subtracted to estimate the elevation change rates Δht for the respective periods. As a time stamp for the SRTM elevation data, the mean date of the SRTM 2000-02-16 is taken. Voids in the SRTM DEM were filled using elevation information from different sources (NASA JPL, 2013) without providing the corresponding date information. Thus, the provided masks used are rejected the void-filling data in our Δht analyses. Since, the glaciers in Bolivia show in general a retreat pattern (Cook and others, Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016), the glacier outlines of the starting year of the respective periods are applied to measure Δht and consequently ice volume losses on the glacierized areas. Additionally, pixels with slopes >50° (3.9% of the glacierized area) are not considered in the further analyses since DEMs are less reliable on steep slopes (Toutin, Reference Toutin2002) and nearly no ice is aggregated there (based on field observations). Gaps in the Δht data on glaciers are filled by hypsometric extrapolation at 100 m elevation bands. In order to remove outliers, for each elevation band, the Δht values are filtered by applying three times the normalized median absolute deviation (NMAD) (Brun and others, Reference Brun, Berthier, Wagnon, Kääb and Treichler2017) and the resulting filtered mean Δht values are applied for the individual elevation bands. The void-filled SRTM DEM is used as the elevation reference for the hypsometric analyses.

Information on the regional glacier mass budget is revealed by calculating the geodetic mass balance (ΔMt). This is defined by Fountain and others (Reference Fountain, Krimmel and Trabant1997) as the elevation change rate integrated over the glacierized region multiplied by an average density ρ, also called volume to mass conversion factor. Two density scenarios were applied. For the first one, an average density of 850 kg m−3 is considered, as suggested by Huss (Reference Huss2013) for alpine glaciers. The second scenario uses two densities of 600 and 900 kg m−3 for areas above and below the ELA (5144 m a.s.l., Rabatel and others, Reference Rabatel2012), respectively (e.g. Gardelle and others, Reference Gardelle, Berthier and Arnaud2012; Kääb and others, Reference Kääb, Berthier, Nuth, Gardelle and Arnaud2012). Kääb and others (Reference Kääb, Berthier, Nuth, Gardelle and Arnaud2012) suggested that the second density scenario is more suitable for glaciers where the volume changes are mainly forced by ablation and accumulation and not by glacier dynamics. Therefore, we used the results of the second scenario for the further analyses and discussion and kept the results of the first scenario for comparison with other studies, which used a similar density approach. Huss (Reference Huss2013) showed that the density can vary significantly. According to these findings, an uncertainty on the density of ±60 kg m−3 is applied for mass-balance quantification over periods longer than 10 years and mass changes larger than ±200 kg m−2 a−1. Otherwise, an uncertainty of ±300 kg m−3 is taken (see also Braun and others, Reference Braun2019), that is, for the period 2013–2016.

The geodetic mass balance of the individual glaciers is computed according to Cogley and others (Reference Cogley, Hock, Rasmussen, Arendt, Bauder, Braithwaite, Jansson, Kaser, Möller, Nicholson and Zemp2011), and the voids in the elevation change fields are handled following the recommended method by McNabb and others (Reference McNabb, Nuth, Kääb and Girod2019). For each glacier with <60% voids in the Δht fields, the hypsometric distribution of Δht is computed. Elevation intervals of 50 m are applied for glaciers with an elevation range >500 m, otherwise 10% of the glacier's elevation range is used for the binning. Information for elevation bands without Δht data is obtained by fitting a third-order polynomial to the hypsometric Δht distribution. Only glaciers with a data coverage of at least two-thirds of the elevation bands are considered in the analysis. The resulting mass budgets of the individual glaciers are analyzed in combination with topographic attributes (area, median elevation, aspect) in order to detect correlations. Mass balances of individual glaciers are computed by applying the second density scenario. Artifact in the Δht fields can bias the mass-balance results of individual glaciers. Therefore, a three times NMAD filter is applied in the hypsometric analysis (see regional analysis above) and finally a 2–98% quantile filter is used to discard possible remaining outliers.

3.3 Uncertainty analysis

The accuracy of the area measurements (δ S) is estimated according to Malz and others (Reference Malz2018). The authors used the accuracy estimation of 3% for outlines of alpine glaciers, mapped by means of LANDSAT imagery from Paul and others (Reference Paul2013), and included a scaling factor for the different perimeter–area ratios.

The uncertainty in the geodetic mass-balance estimates (δ ΔMt) is calculated by:

(1)$$\hskip-5pt \scale 88% { \delta _{\Delta {\rm M}/\Delta {\rm t}} = \sqrt {{\left( {\displaystyle{{\Delta M} \over {\Delta t}}} \right)}^2\left( {{\left( {\displaystyle{{\delta_{\Delta {\rm h}/\Delta {\rm t}}} \over {\Delta h/\Delta t}}} \right)}^2 + {\left( {\displaystyle{{\delta_{\rm S}} \over S}} \right)}^2 + {\left( {\displaystyle{{\delta_\rho} \over \rho}} \right)}^2} \right) + {\left( {\displaystyle{{V_{{\rm pen}}} \over {\Delta t}}{^\ast}\rho} \right)}^2}. } $$

It consists of the error contributions of the elevation change measurements (δ Δh/Δt), the glacier inventory (δ S), the ice density (δ ρ) and an uncertainty due to the differences of the SAR signal penetration (V pent). The accuracy of the averaged Δht measurements (δ Δh/Δt) is affected by the relative vertical accuracy over the survey glacier area (δ Δh/Δt(vert)) and the uncertainty caused by the hypsometric extrapolation. The relative vertical error δ Δh/Δt(vert) is evaluated by analyzing Δht on the ice-free areas and masking out water bodies and vegetation (see Section 3.2). Strongly deviation values are removed by applying a 2–98% quantile filter. Since the offsets depend on the slope (Fig. 2; Gardelle and others, Reference Gardelle, Berthier and Arnaud2012; Vijay and Braun, Reference Vijay and Braun2016), the area weighted standard deviations (σ AW) based on the slope distribution (5° bins) are calculated. Additionally, the deviations for each slope bin are filtered by applying 3*NMAD before computing σ AW. To estimate the accuracy of the glacier-wide averaged elevation change (δ Δh/Δt(vert)), we followed the approach of Rolstad and others (Reference Rolstad, Haug and Denby2009) (applied by, e.g. Fischer and others, Reference Fischer, Huss and Hoelzle2015; Brun and others, Reference Brun, Berthier, Wagnon, Kääb and Treichler2017) in order to take into account spatial out correlation. For glacier areas (A gl) larger than the correlation area (A cor = π × d2cor; dcor: correlation length) δ Δh/Δt(vert) is calculated as:

(2)$$\delta _{{\rm \Delta h/\Delta t}\lpar {{\rm vert}} \rpar } = \sqrt {\displaystyle{{A_{{\rm cor}}} \over {5A_{{\rm gl}}}}} \sigma _{{\rm AW}}.$$

Fig. 2. Area on (light blue) and off glacier (red) and off glacier dh/dt distribution (blue dots) based on slope. Error bars indicate NMAD of dh/dt values for each slope bin. Dotted line shows the slope threshold used in the analysis. Note: the glacier area is scaled by a factor of 10.

For A gl<A cor as:

(3)$$\delta _{{\rm \Delta h/\Delta t}\lpar {{\rm vert}} \rpar } = \sigma _{{\rm AW}}.$$

Since we are integrating the elevation change over large regions for the regional analysis, the above approach is applied for each contiguous glacierized area and not for individual glaciers. For the glacier-specific mass balances, the individual glacier areas are used. The region-wide δ Δh/Δt(vert) is calculated as the area weighted average of all individual glacierized areas. The spatial auto-correlation of our off-glacier elevation change measurements is analyzed by generating semi-variograms with 1 00 000 random samples, binned in 30 m distance intervals for a maximum distance of 20 km. The correlation range is estimated by applying a spherical semi-variogram function and a mean d cor of 394 m is found for our data sets (an example is illustrated in Fig. S1).

In order to account for the error due to the hypsometric extrapolation, we followed the approach of Berthier and others (Reference Berthier2014) and multiplied the error δ Δh/Δt(vert) by a constant factor (we chose a factor of 2, Braun and others, Reference Braun2019) for the unmapped areas to finally retrieve δ Δh/Δt.

SAR signal penetration is evaluated by comparing a DEM derived from a Pléiades acquisition on 28 May 2013 (see Fig. 1) and the TDX DEM in 2013 (TDX data covering the Pléiades DEM is from 31 January 2013). Both DEMs are coregistered (see Section 3.2) and the hypsometric distribution of the elevation differences in glacier areas (inventory from 2013) are analyzed (Fig. S2). A negligible influence on the elevation change measurements by the SAR signal penetration is found (see Section 5.2 for a detailed discussion). Consequently, no correction of the SAR data is applied. However, we considered a potential influence in the total error quantification of the mass balances, since the surface conditions in 2000 and 2016 might have been different as in 2013 and SRTM data were acquired at a lower SAR frequency.

The uncertainty contribution due to an offset in the SAR signal penetration between SRTM and TanDEM-X is estimated according to Malz and others (Reference Malz2018). Deviations in the elevation change rates caused by differences in the SAR signal penetration of C- and X-band data in snow and ice surfaces in areas below the ELA (estimated at an average elevation of 5144 m a.s.l. for Zongo Glacier by Rabatel and others (Reference Rabatel2012), and assumed to be representative of the regional average in the current study) are neglected, since no penetration is assumed in the typical snow and ice melt affected surfaces. Above the ELA, a linear increase toward 5 m of penetration at the maximum elevation of the study region is applied to estimate V pen. The respective density of the applied scenario is used to convert V pen into mass. Moreover, for the second scenario, the potential offset in the mass budgets caused by differences in the SAR signal penetration is ~40% lower due to the lower density considered above the ELA.

4. Results

4.1 Area changes

The glacier outlines obtained for 2000, 2013 and 2016 are plotted in Figure 3 (detailed subsets are plotted in Fig. S3). A reduction of the extent of the glacierized areas is clearly visible. This shrunk from 281 ± 14 km2 in 2000 to 246 ± 14 km2 in 2013 (−12%) and 200 ± 11 km2 in 2016 (−29%) amounting to an area loss of 81 ± 18 km2. The area change values of the study area are summarized in Table 1. The average shrinkage rates show a stronger recession toward recent years from −2.7 ± 1.5 to −15 ± 5 km2 a−1 in the periods 2000–2013 and 2013–2016, respectively. Moreover, the number of glaciers >0.01 km2 reduced from 332 in 2000 to 272 glaciers in 2016 and 34 glaciers (nine glaciers in 2000–2013) completely vanished. Additionally, glacier fragmentation was observed at several glaciers (see Fig. S3). In the Cordillera Real, the glacier area shrunk from 244 km2 in 2000 to 175 km2 in 2016 (−28%) and in Tres Cruces from 37 km2 in 2000 to 24 km2 in 2016 (−34%). The links between the area changes of the individual glaciers (in %) and topographic variables, such as area, median elevation and aspect, are illustrated in Figure 4 for the period 2000–2013 (Figs S4 and S5 for the periods 2000–2016 and 2013–2016). Highest retreat rates are generally found for low-lying and small glaciers. Area change measurements for small glaciers are more sensitive to mapping errors due to the typically higher perimeter–area ratio. This partially leads to a wider range and even positive values in relative area changes for smaller glaciers as compared to the larger ones.

Fig. 3. Glacier area changes. Ice divides (black polygons) are from the glacier outlines in 2000. Dashed polygons indicate subsets illustrated in Figure S3 (Supplement) Background: SRTM hillshade © NASA 2000.

Fig. 4. Relative area changes (dS) (2000–2013) of individual glaciers (dot color) plotted against glacier size (dot size), median elevation (distance from center) and mean aspect (orientation). Red circle: equilibrium line altitude (ELA) from Rabatel and others (Reference Rabatel2012).

Table 1. Observed glacier changes for different time intervals

Δt: observation period; ΔS: glacier area change; ΔSt: glacier area change rate; Δh Mt: average measured surface lowering rate on S M; Δh Et: average extrapolated surface lowering rate on S; S: analyzed glacier area at the beginning of the observation period. Note: slopes >50° are excluded; S M: fraction of S covered by Δh Mt measurements; ΔMt: mass balance (average and specific); ΔM: total mass change in the observation period.

a Volume to mass conversion factor of 850 kg m−3.

b Volume to mass conversion factor of 900 and 600 kg m−3 for areas below and above ELA, respectively.

4.2 Elevation changes and mass balance

Figure 5 shows the measured (unfiltered) surface elevation change rates for the periods 2000–2013 and 2013–2016 (detailed subsets are presented in Figs S6 and S7). A clear surface lowering at most of the glacier termini is found. Gaps in the Δht fields are due to incomplete SRTM data coverage (coverage of 67% of the glacierized area) and voids in the TDX DEMs (coverage of 94 and 92% of the glacierized areas in 2013 and 2016, respectively, within the TDX footprints) caused by typical SAR limitations like layover and shadowing (see also Sections 3.2 and 5.2). In the period 2000–2016, the elevation change is measured only on 50% of the glacier area since no complete coverage of the study site by TDX data is available in 2016 (Fig. 1). The hypsometric distribution of the glacier area and the obtained Δht (filtered) are illustrated in Figure 6 for the period 2000–2013 (Figs S8 and S9 for the periods 2000–2016 and 2013–2016, respectively). Most of the surface lowering is observed at elevations below the ELA. The measured surface elevation gains in the upper glacier areas can be neglected considering the small fraction of coverage (<1%). An increase in the extrapolated mean surface elevation change from −0.51 ± 0.11 m a−1 for 2000–2013 to −0.74 ± 0.34 m a−1 for 2013–2016 is found. All mean Δht estimates are presented in Table 1. Consequently, a more negative geodetic mass balance of −0.12 ± 0.09 Gt a−1 is revealed for the period 2013–2016 as compared to −0.10 ± 0.03 Gt a−1 in the period 2000–2013. In total, −1.8 ± 0.5 Gt of ice was lost in the observation period.

Fig. 5. Glacier surface elevation changes (unfiltered; left panel: 2000–2013; right panel: 2013–2016). Dashed polygons indicate subsets illustrated in Figures S6 and S7 (Supplement). Background SRTM hillshade © NASA 2000.

Fig. 6. Hypsometric distribution of glacier area (light blue), glacier area with dh/dt measurements (red) and mean, filtered dh/dt values (blue dots) of each hypsometric bin for the observation period 2000–2013. Error bars indicate NMAD of dh/dt for each hypsometric bin.

Glacier-wide mass balances for ~60% of all glaciers could be quantified for the study periods 2000–2013 and 2013–2016. Due to voids in the SRTM DEM and incomplete coverage by the TDX DEM in 2016, individual glacier mass balances of only 40% of all glaciers are obtained in the period 2000–2016. Histograms indicating the distribution of the glacier-specific mass balances, as well as minimum and maximum values, are provided in Figure S10. Figure 7 illustrates links between the mass balances obtained for each glacier and topographic variables for the period 2000–2013 (Figs S11 and S12 for the periods 2000–2016 and 2013–2016, respectively). A trend of more negative mass balances for glaciers with lower median elevation as well as for glaciers with northward orientation is revealed.

Fig. 7. Specific mass balance (spMB) (2000–2013) of individual glaciers (dot color) plotted against glacier size (dot size), median elevation (distance from center) and mean aspect (orientation). Red circle: equilibrium line altitude (ELA) from Rabatel and others (Reference Rabatel2012). Note: only glaciers with >40% elevation change data coverage, which is spread over >2/3 of the hypsometric distribution are included.

5. Discussion

5.1 Area changes

Our quantified glacier area loss of up to 30% (2000–2016) in the study region is in line with previous observations in this region and highlights the dramatic glacier changes in Bolivia. Several studies (e.g. Rabatel and others, Reference Rabatel, Machaca, Francou and Vincent2006; Soruco and others, Reference Soruco2009b; Réveillet and others, Reference Réveillet, Rabatel, Gillet-Chaulet and Soruco2015) have reported glacier recession of individual glaciers in the Cordillera Real. Soruco and others (Reference Soruco, Vincent, Francou and Gonzalez2009a) reported in the central part of the Cordillera Real an area change of −48% (−1.5% a−1) between 1975 and 2006 and −1.7% a−1 between 1997 and 2006, which is comparable to our observed retreat rate in the Cordillera Real of −28% (−1.8% a−1) in the period 2000–2016. In Tres Cruces, Albert and others (Reference Albert, Kargel, Leonard, Bishop, Kääb and Raup2014), Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) and Veettil and others (Reference Veettil, Wang, Simões and Pereira2018) measured glacier area changes of ~−55% (−1.6% a−1) for 1975–2009, −47.3% (−1.7% a−1) for 1986–2014 and −35.1% (−1.7% a−1) for 1995–2016, respectively. We observed a higher retreat rate of −2.1% a−1 for the period 2000–2016 in Tres Cruces, which agrees with the reported increased glacier retreat in Tres Cruces after 2010 by Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016). Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) mapped area changes throughout the Bolivian Cordillera Oriental for the first time and reported recessions of 43.1% (1.5% a−1) for the period 1986–2014, which is similar to the observation of 1.4% a−1 by Veettil and others (Reference Veettil, Wang, Simões and Pereira2018) for 1995–2016. Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) found high retreat rates of 8.1 km2 a−1 (2.5% a−1) in the period 1999–2010 and lower rates of 4.0 km2 a−1 (1.2% a−1) in the period 2010–2014, including the Cordillera Apolobamba. Taking into account the different sizes of the studied areas (351 km2 in 2014, including the Cordillera Apolobamba), this leads to an average loss rate of 1.4% a−1 for 1999–2014, based on the results of Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016). This value is comparable to our observed retreat rates of 0.9% a−1 and 1.8% a−1 in the periods 2000–2013 and 2000–2016, respectively, considering the shift in the observation intervals and a bias by the higher retreat rates in the Cordillera Apolobamba as compared to the Cordillera Real (observed by Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) and Veettil and others (Reference Veettil, Wang, Simões and Pereira2018)).

In the time interval 2013–2016, we obtained an enhancement of the retreat rates toward 15 ± 5 km2 a−1 in our study area. A similar area change trend was observed at Peruvian glaciers by Seehaus and others (Reference Seehaus2019). The increased rates are comparable to the finding of Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) in the period 1986–1992 (14.5 km2 a−1). These high retreat rates are most likely to be associated with the strong El Niño events that occurred during these two observation periods (see Section 5.3 and Seehaus and others (Reference Seehaus2019) for more details).

Our mapped glacier area for 2000 is similar to that measured by Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016). It is ~13.3% smaller than the glacier extents based on the RGI6.0 (2000–2003), which was expected according to the limitation indicated in the Technical Report of the RGI 6.0 (see Section 3.1). Moreover, the minimum glacier altitude of 3706 m a.s.l. according to the RGI 6.0 is caused by a misclassified glacier area. We find a more reasonable minimum glacier altitude of 4418 m a.s.l. Figures 4, S4 and S5 indicate that the highest relative area changes are in general found for small and low-lying glaciers, confirming the findings by Veettil and others (Reference Veettil, Wang, Simões and Pereira2018). This can be partly attributed to the fact that small glaciers are generally located at lower elevations and thus experience more ablation due to changes in the climatic conditions (see Section 1). Some of the small and low-lying glaciers even lost their accumulation areas due to the shift of the ELA (Vuille and others, Reference Vuille2008). Rabatel and others (Reference Rabatel2013) reported the disappearance of the monitored Chacaltaya Glacier in the Cordillera Real in 2009. Its extinction was already projected by Ramirez and others (Reference Ramirez2001) since the estimated mean ELA was above its upper extent. The authors proposed that Chacaltaya Glacier is representative of the Bolivian glaciers, of which ~80% are smaller than 0.5 km2 (based on our glacier inventories, 63 and 76% were smaller than 0.5 km2 in 2000 and 2016, respectively) with average altitudes close to or below the ELA. All of our identified extinct glaciers (2000–2016) had a median elevation below 5400 m a.s.l. and were smaller than 0.5 km2, which supports the predictions of Ramirez and others (Reference Ramirez2001). At Cerro Charquini (5392 m a.s.l.), Rabatel and others (Reference Rabatel, Machaca, Francou and Vincent2006) derived an upward shift of the ELA by 160 m between the LIA maximum (~1660 AD) and 1997 and attributed it to a temperature increase by ~0.6 °C as the main driving factor. For Zongo Glacier, Rabatel and others (Reference Rabatel2012) estimated an average ELA of 5144 ± 67 m a.s.l. for the period 1991–2006 and Vuille and others (Reference Vuille2018) computed an ELA shift toward 5400 and 5700 m a.s.l. by 2100 based on the CIMP5 RCP4.5 and RCP8.5 scenarios, respectively. Consequently, further disappearance of glaciers can be expected considering the large number of glaciers at low altitudes with area losses of more than 60% (2000–2016) and the projected upward shift of the ELA.

The higher shrinkage rates at Tres Cruces, as also observed by Veettil and others (Reference Veettil, Wang, Simões and Pereira2018), can be related to the observed link between median glacier elevations and shrinkage as well. The median glacier elevation in the study region for the 2000 inventory was 5342 m a.s.l. (5365 m a.s.l. for Cordillera Real only), whereas at Tres Cruces the median glacier elevation was 5285 m a.s.l. Moreover, the uppermost elevation at Tres Cruces is 5728 m a.s.l., clearly below the maximum elevation at Cordillera Real of 6425 m a.s.l. Thus, the glaciers at Tres Cruces are even more strongly affected by the changing climate conditions and the associated uplift of the ELA than the glaciers in Cordillera Real.

Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) have evaluated the formation of proglacial lakes and the glacier lake outburst flood (GLOF) risks in Bolivia from 1986 to 2014. The continuous large-scale glacier shrinkage revealed in this region until 2016 might have led to further extension or formation of new proglacial lakes. Moreover, further rapid glacier shrinkage is expected (e.g. Ramirez and others, Reference Ramirez2001; Vuille and others, Reference Vuille2018). Thus, we strongly support the suggestion by Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) to further monitor proglacial lakes and GLOF risks in Bolivia.

5.2 Elevation changes and mass balances

An average glacier surface elevation change of −9.9 m (−0.60 m a−1) and highest surface lowering of up to ~50 m (~3.1 m a−1) at the lower glacier reaches are derived from SAR remote-sensing data for the study period 2000–2016. The estimated accuracy of the elevation change (δ Δh/Δt) quantified over stable terrain of 0.88–1.04 m (long-term average: 0.08 m a−1 for 2000–2016; short-term average: 0.24 m a−1 for 2013–2016) is comparable to or even better than the reported values of similar studies in other mountain regions (e.g. 0.22–0.28 m a−1 for 2000–2012 by Vijay and Braun (Reference Vijay and Braun2018), 0.43 m a−1 for 2000–2013 by Vijay and Braun (Reference Vijay and Braun2016), 0.22 m a−1 for 1999–2008 by Gardelle and others (Reference Gardelle, Berthier and Arnaud2012) and 0.21–0.28 m a−1 for 2000–2015/16 by Malz and others (Reference Malz2018)). The evaluation of the elevation differences on glacier areas between the TDX DEM and Pléiades DEM in 2013 revealed an average offset at regions below and above the ELA of −0.7 and −0.1 m, respectively (Fig. S2). The total offset due to SAR signal penetration is small considering the temporal difference of 4 months between both acquisitions (end of wet season and beginning of dry season, accumulation in the upper reaches and high melt rates toward the terminus), the average surface-lowering rates (2000–2016) of −1.0 and −0.7 m a−1 at areas below and above the ELA, respectively, and the uncertainty of the elevation change measurements. The SAR data (SRTM and TDX) were acquired during the transition and wet season (Sections 2 and 3.2), when the melt rates are highest and meltwater is most likely present in the snow pack, strongly limiting the signal penetrations. Moreover, we compared only DEMs based on SAR acquisitions, thus only relative differences in the SAR signal penetration affect the volume change estimation. Therefore, we conclude, that the impact of the SAR signal penetration on elevation change computations is negligible. However, in order to be consistent with previous studies and to account for potential differences in the glacier surface conditions in 2000 and 2016 as well as for the differences in the SAR frequency, a contribution to the error budget of the mass balances is considered. The estimated uncertainty V pen, due to the difference in the SAR signal penetration between SRTM and TanDEM-X (C-band vs X-band), amounts to 0.288 km3 (area above ELA: 231 km2, ~85.5% of the surveyed area) for the periods 2000–2013 and 2000–2016 and 0.273 km3 (area above ELA: 210 km2, ~88.6%) for the period 2013–2016. The corresponding mean area weighted penetration depth uncertainty is equal to 0.065 m a−1 for the period 2000–2016, which is comparable to the uncertainty of 0.081 m a−1 reported by Malz and others (Reference Malz2018) for a similar observation period (15.86 years).

A strongly negative mass balance of −1030 ± 830 kg m−2 a−1 for the RGI6.0 region ‘Low Latitudes’ (Tropics) is reported by Zemp and others (Reference Zemp2019) for the period 2006–2016. However, this estimate is based on a very sparse data coverage (<2%), leading to the large uncertainty. The authors expressed the need for more observations in the Tropical Andes, which will be consequently provided by this study. Braun and others (Reference Braun2019) observed a mass balance of −340 ± 80 kg m−2 a1 (ρ = 850 kg m−3) for the southern wet outer tropics in 2000–2013, which corresponds to ~22% lower mass loss rate than revealed in this study (−437 ± 122 kg m−2 a−1, first density scenario). As discussed in Section 5.1, the RGI6.0 glacier inventory, which was used by Braun and others (Reference Braun2019), contains many misclassified glacier areas that consequently led to lower specific mass loss. The analysis at the central Cordillera Real by Soruco and others (Reference Soruco, Vincent, Francou and Gonzalez2009a) revealed a volume loss of 43% in the period 1963–2006 based on the elevation change measurements at 21 glaciers. The authors observed also large differences for the measured glaciers, ranging from −260 to −1380 kg m−2 a−1 (ρ = 900 kg m−3). Moreover, they found strong temporal variations ranging from an average mass balance of 700 kg m−2 a−1 for 1963–1975 to −500 kg m−2 a−1 for 1975–1983. For the period 1997–2006, an average mass balance of −450 kg m−2 a−1 was derived for the surveyed glaciers. This value is similar to our observations of −437 ± 122 and −499 ± 97 kg m−2 a−1 for the periods 2000–2013 and 2000–2016, respectively (first density scenario). However, for the period 2013–2016, we measured a higher mass loss of −630 ± 450 kg m−2 a−1, which correlates well with our observed higher glacier retreat rates for this time interval (Fig. 8). An average glaciological mass balances at Zongo Glacier of −279 and −273 kg m−2 a−1 is reported for the period 2000–2013 and 2000–2016, respectively (WGMS). These values are comparable to our findings of −275 ± 70 and −365 ± 60 kg m−2 a−1 in the respective intervals at Zongo Glacier. At Charquini Sur Glacier, we quantified a mass balance of −839 ± 160 kg m−2 a−1 in the period 2000–2016, which fits the average measured glaciological mass balance of −967 kg m−2 a−1 in the period 2002–2016 (WGMS). Chacaltaya Glacier disappeared within the observation period. Therefore, a comparison with glaciological mass-balance estimates is not meaningful.

Fig. 8. Temporal evolution of glacier area and specific mass balance (spMB) (a and b) as well as Oceanic Nino Index (ONI, c) for the studied period.

At Zongo Glacier, Soruco and others (Reference Soruco2009b) found that 80% of the specific mass-balance contribution comes from the strongly ablating low-elevation area (about one-third of the glacier area). This fits well with our observed pattern on elevation change and altitude (Figs 6, S8 and S9) and also explains the high retreat rates of the glacier tongues. Figures 7, S11 and S12 indicate that small, low-lying and northward oriented glaciers showed predominantly more negative mass balances. Similar to Soruco and others (Reference Soruco2009b), we correlated the glacier orientation and mass budgets. The applied sinusoidal fits show minimum mass balances for North/North-West orientations (Fig. S14), which is in accordance with observations by Soruco and others (Reference Soruco2009b) who reported the highest mass losses on glaciers facing toward North/North-West since the annual mean incoming solar radiation is highest for north-facing slopes at these latitudes. In particular, some of the smallest glaciers show highest specific mass loss variability since they are strongly affected by local settings (Vuille and others, Reference Vuille2008). The air temperature at surrounding rock surfaces can rise above 20 °C and lead to advection of warm air over the glacier (Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003). This so-called ‘edge effect’ becomes important once the glacier shrunk below a critical size and can foster accelerated glacier recession. This pattern also fits with our observed higher relative area changes for small low-lying glaciers (Section 5.1).

Rabatel and others (Reference Rabatel2013) derived two times higher mass loss rates for glaciers with maximum elevations below 5400 m a.s.l. as compared to glaciers with maximum altitudes above this limit, based on mass-balance series of eight glaciers throughout the Tropical Andes (three glaciers in Bolivia). Such a distinct pattern is not so obvious in the geodetic measurements obtained in this study. However, we revealed an overall trend toward higher mass loss rates for glacier with a lower upper limit, and 12–24% stronger mass losses for glaciers with maximum altitudes below 5400 m a.s.l. (Fig. S13). On the one hand, our analysis covers only Cordillera Real and Tres Cruces, on the other hand, the representativeness of the eight glaciers, analyzed by Rabatel and others (Reference Rabatel2013), for all glaciers in the Tropical Andes is not guaranteed. Moreover, glaciological mass-balance measurements can be biased by the spatial data sampling (see above). Thus, we attribute the offset in the trends to differences in the spatial coverage and sampling as well as to the different mass-balance method used.

As mentioned in Section 4.2, data voids in InSAR-based DEMs can occur due to topographic limitations like SAR shadow and layover. For the study periods 2000–2013 and 2000–2016, the SRTM DEM is used as a reference. Only 63.3% of the glacierized area in 2000 is covered by SRTM elevation data. In combination with the TanDEM-X coverage in 2013, we retrieve Δht measurements on 59.2% of the glacier area. This corresponds to coverage of 93.5% of the SRTM data by TanDEM-X data. The largest data voids are present in the Δht fields for 2000–2016 (52.0%) since no complete TanDEM-X coverage of our study site is available in 2016 (~47 km2 ≙ 20% of the glacier area in 2013 is not covered, Fig. 1). However, the amount of voids in the Δht data does not significantly influence the region-wide mass-balance estimates when applying ‘global hypsometric interpolation’. According to McNabb and others (Reference McNabb, Nuth, Kääb and Girod2019), this method is suitable to compute regional averages from Δht data sets with up to 60% of voids. For the study period 2013–2016, the mosaicked TanDEM-X DEM from 2013 is used as a reference, leading to the highest coverage of 69.2% (62.3% after filtering) of all analyzed glacierized areas by Δht data. Moreover, coverage of 86.8% by Δht data on the glacier areas within the TanDEM-X footprints in 2013 and 2016 (regions where TanDEM-X data from both years is available, see also Fig. 1) is revealed for 2013–2016. Consequently, the achievable data coverage can be increased when comparing TanDEM-X to TanDEM-X DEMs, especially in alpine regions where SRTM data can have large voids. Future studies can benefit from using TanDEM-X DEMs instead of SRTM DEM as a reference since less extrapolation is needed, helping to increase the accuracy (Section 3.3). Moreover, the amount of individual glacier mass-balance estimates is limited by voids in the Δht fields (Section 3.2). Thus, a higher spatial coverage will also increase the amount and quality of glacier-wide mass budget information. Combining TanDEM-X data from ascending and descending orbits is highly recommended if available for the specific region and time in order to obtain optimum data coverage. Moreover, many published studies used SRTM DEMs to compute elevation change rates. However, it is often not clear how the void-filled areas (no timestamp is available) were handled. Therefore, we emphasize to clearly document for each study which SRTM DEM version is used and especially how the void-filled areas are handled when using SRTM DEMs to assess elevation change rates.

5.3 Climatic effects on glacier changes

Previous surveys at Zongo and Chacaltaya glaciers (e.g. Wagnon and others, Reference Wagnon, Ribstein, Francou and Sicart2001; Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003) demonstrated that annual surface mass balance is dominated by ablation and accumulation processes during the wet summer season (70% of the total variance at Zongo Glacier). For the observation period 2000–2013, the glacier mass balance is computed by comparing DEMs from the same seasons (around February, see Table 1) to avoid a potential bias due to intra-annual mass-balance fluctuations and to minimize the impact of differences in the SAR signal penetration due to different surface conditions. Unfortunately, in 2016, TanDEM-X data are only available for the period September–November (see Table 1), the end of the dry season and the beginning of the wet season. Minimal precipitation occurs typically during the dry season and the glacier mass balance is dominated by ablation (Kaser, Reference Kaser2001; Favier and others, Reference Favier, Wagnon and Ribstein2004). Thus, we would like to point out that the presented mass balances for 2000–2016 and 2013–2016 could potentially be slightly more positive when accounting for the seasonal bias. However, it is difficult to correctly quantify such a bias. Therefore, no seasonal correction is applied.

The inter-annual variability of the surface mass balance is strongly influenced by ENSO events on regional scales (Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003) (see also Section 1). In order to evaluate the ENSO effects on the glacier changes, we used as an indicator the monthly Oceanic Niño Index (ONI) provided by NOAA Climate Prediction Center (http://origin.cpc.ncep.noaa.gov/products/analysis_monitoring/ensostuff/ONI_v5.php). According to NOAA's definitions, El Niño is present for ONI values above +0.5 and La Niña conditions are indicated for ONI values below −0.5. The temporal evolution of ONI is plotted in Figure 8 for the period 2000–2016. A strong El Niño event is clearly visible for 2015 (ONI values of up to 2.6), which leads to a clear positive average ONI value of 0.42 in the period 2013–2016. In the period 2000–2013, the ONI fluctuations are more balanced with a mean ONI value of −0.17. More negative mass balances, caused by precipitation deficits and higher incoming radiation energy, are typical for El Niño events in the Bolivian Andes (Francou and others, Reference Francou, Vuille, Wagnon, Mendoza and Sicart2003; Maussion and others, Reference Maussion, Gurgiser, Großhauser, Kaser and Marzeion2015). A two-thirds higher runoff at Zongo Glacier was observed by Wagnon and others (Reference Wagnon, Ribstein, Francou and Sicart2001) during the strong El Niño event in 1997–1998. Additionally, Soruco and others (Reference Soruco, Vincent, Francou and Gonzalez2009a) revealed that 60% of the mass loss in 1997–2006 occurred during this event. Schoolmeester and others (Reference Schoolmeester2018) reported an extreme drought that stressed the water availability in Bolivia in 2015/16. Moreover, we analyzed monthly mean temperature recordings at Zongo and Charquini Sur glaciers (Figs S15 and S16). Higher average temperatures are obtained for the period 2013–2016, as compared to the period 2000–2013, and a clear peak in the temperature trend is obvious during the El Niño event in 2015/16. Thus, it is likely that our observed increase glacier shrinkage in the period 2013–2016, which is also reported by Seehaus and others (Reference Seehaus2019) for Peruvian glaciers, can be attributed to the strong El Niño event in 2015/16. Since the response time of tropical glaciers to mass-balance variations is very short (e.g. Francou and others, Reference Francou, Vuille, Favier and Cáceres2004; Rabatel and others, Reference Rabatel2013) and the fact that the glaciers in the low latitudes are in average the thinnest worldwide (Farinotti and others, Reference Farinotti2019), the increased melt in the period 2013–2016 explains the enhanced glacier recession in this interval. In the period 1986–1992, a mean ONI value of +0.24 (maximum values of 1.71 in January 1992 and 1.70 in August 1987) was obtained. Consequently, the very high glacier retreat rates obtained by Cook and others (Reference Cook, Kougkoulos, Edwards, Dortch and Hoffmann2016) of 14.5 km2 a−1 in this period, which are similar to ours in the El Niño period 2013–2016, can be most likely attributed to increased El Niño activities. This link is in line with the observations by Morizawa and others (Reference Morizawa, Asaoka, Kazama and Gunawardhana2013) who derived higher recession rates during El Niño and even partially slight area gains during La Niña at Condoriri Glacier, Cordillera Real in the period 1988–2010. Moreover, Silverio and Jaquet (Reference Silverio and Jaquet2017) revealed a similar pattern and presented a linear correlation (R 2 = 0.8) between ONI and glacier area changes. We correlated glaciological measurements at Zongo and Charquini Sur Glacier (Figs S17 and S18, WGMS) with average ONI values of the respective analysis period (September–August). Tendencies toward smaller accumulation-area-ratios (AAR), higher ELA values and more negative specific mass balances are found (Figs S19 and S20). These revealed tendencies fit the observations by the previously mentioned studies and support our assumption that the strong El Niño event in 2015/16 lead to the increased glacier wastage.

Vuille and others (Reference Vuille2018) estimated an average temperature increase by 1°–5 °C until 2100 relative to the 1961–1960 mean value for the Tropical Andes, based on CMIP5 scenarios (RCP 4.5 and RCP 8.5). The authors also evaluated the impact on the ELA and projected for Zongo Glacier a rise toward ~5400 m a.s.l. (RCP4.5) and ~5700 m a.s.l. (RCP8.5). Taking into account our observed median glacier elevations and the relation of area changes to median glacier elevation (Section 5.1), such warming will further accelerate the glacier recession and lead to the disappearance of many small and low-lying glaciers and lead to significant volume losses as projected by Réveillet and others (Reference Réveillet, Rabatel, Gillet-Chaulet and Soruco2015), for example, 40–89% at Zongo Glacier (depending on the model scenario). The impact of future glacier recession on the basin runoff was evaluated by Huss and Hock (Reference Huss and Hock2018) on global scales. Their estimates indicated that the Bolivian glaciers have already reached or are close to ‘peak water’ (the point in time with maximum annual glacier runoff) and that a reduction of future glacier runoff in the range of 20–40% until 2090 has to be expected for the major drainage basins of the Bolivian glaciers, the Amazon and Titicaca basins. Comparable reduction in glacier meltwater contribution to the runoff of Zongo watershed by 2100 was modeled by Frans and others (Reference Frans2015). Large-scale and spatially detailed information is needed on the future evolution of glaciers, since the metropolis La Paz/El Alto obtains water from a large catchment area and further extensions are planned (personal communication with local stake holders). Thus, the results obtained in this study provided a consistent database for further regional and detailed analyses of glacier changes in Bolivia.

6. Conclusions

The results presented in this study clearly indicate the rapid glacier shrinkage in Bolivia. We observed area changes of −81 ± 18 km2 (−29%, −5.1 km a−1), the disappearance of 34 glaciers (reduction to an area <0.01 km2) and an average mass-balance rate of −403 ± 78 kg m−2 a−1 (−1.8 ± 0.5 Gt) in the period 2000–2016. An increased area change rate of −15 ± 5 km2 a−1 is revealed for 2013–2016 and can be attributed to the strong El Niño event in this period. The more negative mass-balance rate of −488 ± 349 kg m−2 a−1 in this period (not significant, considering the uncertainties) and reported values for Peruvian glaciers support these findings. The revealed temporal variability of glacier changes highlights the impact of short-term climatic variations on glacier mass balances in the Tropical Andes. The analysis of the area and mass changes of individual glaciers revealed a trend of higher mass losses for northward oriented glaciers and increased glacier wastage for small glaciers at low altitudes confirming the predicted vanishing of small low-lying glaciers in this region. Our observations also confirm that glacier retreat rates are even higher in the Cordillera Tres Cruces than in the Cordillera Real, which we attributed to the lower average glacier elevations.

Additionally, we showed that gaps in the SRTM are the major source of voids in the glacier surface elevation change fields and that by using solely TanDEM-X data, a coverage of 78% could be achieved.

This study provides the first multi-temporal, region-wide and spatially detailed mass-balance estimation for Cordillera Real and Cordillera Tres Cruces combined with a temporally consistent analysis of glacier area changes. These findings form fundamental information for the modeling of future glacier evolution, for developing water resource management plans and for further monitoring of glaciers. The ongoing dramatic glacier wastage in Bolivia will likely cause serious socio-economic challenges in the next decades. On the one hand, it will lead to a reduction in glacier meltwater contribution to discharge and hence affect regional water availability for mining, hydropower production, irrigation or domestic water supply. On the other hand, it can foster the extension and formation of glacial lakes, leading to a potential increase in the risk of GLOFs endangering downstream communities. Thus, further detailed monitoring of the Bolivian glaciers is highly advisable from scientific but also due to socio-economic aspects.

Supplementary material

The supplementary material for this article can be found at https://doi.org/10.1017/jog.2019.94

Acknowledgements

This work was financially supported with the DLR/BMWi grant GEKKO (50EE1544), the grant BR2105/14-1 within the DFG Priority Program ‘Regional Sea Level Change and Society’ as well as the HGF Alliance Remote Sensing & Earth System Dynamics. We thank the German Aerospace Center for providing TanDEM-X data free of charge under AO XTI_GLAC0264. Landsat data were kindly provided via USGS Earth Explorer; SRTM by NASA. In situ glaciological and meteorological measurements at Zongo and Charquini Sur glaciers were provided through the French Service National d'Observation GLACIOCLIM (https://glacioclim.osug.fr/) and the Laboratoire Mixte International GREAT-ICE (a joint initiative between the French IRD and national counterpart in Ecuador, Bolivia, Peru and Colombia). The Pléiades DEM was obtained through the CNES/SPOT-Image ISIS program (contract #FC18473, led by A. Rabatel). A. Rabatel acknowledges the support of the Labex (Investissements d'avenir – ANR10 LABX56). We thank the scientific editor Hamish Pritchard and four anonymous reviewers.

Author contributions

TS led the study, analyzed the data and wrote the manuscript. TS, PM and CS developed jointly and wrote the analysis routines. AS and AR contributed to the interpretation of the data and provided field measurements and the 2013 Pléiades DEM. MB initiated and supervised the project. All authors revised the manuscript.

Conflict of interest

The authors declare no conflicts of interest. The founding sponsors had no role in the design of the study; in the collection, analyses or interpretation of data; in the writing of the manuscript; and in the decision to publish the results.

Data and materials availability

Elevation change fields are available via the World Data Center PANGAEA operated by AWI Bremerhaven (https://doi.pangaea.de/10.1594/PANGAEA.907325, doi:10.1594/PANGAEA.907325). Glacier area information and glacier-specific results are available via the World Glacier Monitoring Service (WGMS) and GLIMS.

References

Albert, T (2014) In Kargel, JS, Leonard, J, Bishop, MP, Kääb, A and Raup, BH (eds), Global Land Ice Measurements from Space. Springer, Berlin, Heidelberg, 609638 (doi: 10.1007/978-3-540-79818-7_26).CrossRefGoogle Scholar
Berthier, E and 10 others (2014) Glacier topography and elevation changes derived from Pléiades sub-meter stereo images. The Cryosphere 8(6), 22752291. doi: 10.5194/tc-8-2275-2014.CrossRefGoogle Scholar
Braun, MH and 8 others (2019) Constraining glacier elevation and mass changes in South America. Nature Climate Change 1. doi: 10.1038/s41558-018-0375-7.Google Scholar
Brun, F, Berthier, E, Wagnon, P, Kääb, A and Treichler, D (2017) A spatially resolved estimate of High Mountain Asia glacier mass balances from 2000 to 2016. Nature Geoscience 10(9), 668. doi: 10.1038/ngeo2999.CrossRefGoogle Scholar
Buytaert, W and 7 others (2017). Glacial melt content of water use in the tropical Andes. Environmental Research Letters 12(11), 114014. doi: 10.1088/1748-9326/aa926c.CrossRefGoogle Scholar
Cogley, JG, Hock, R, Rasmussen, LA, Arendt, AA, Bauder, A, Braithwaite, RJ, Jansson, P, Kaser, G, Möller, M, Nicholson, L and Zemp, M (2011) Glossary of Glacier Mass Balance and Related Terms. IHP-VII Technical Documents in Hydrology No. IACS Contribution 48(04). doi: 10.1017/S0032247411000805.Google Scholar
Cook, SJ, Kougkoulos, I, Edwards, LA, Dortch, J and Hoffmann, D (2016) Glacier change and glacial lake outburst flood risk in the Bolivian Andes. The Cryosphere 10(5), 23992413. doi: 10.5194/tc-10-2399-2016.CrossRefGoogle Scholar
Farinotti, D and 6 others (2019). A consensus estimate for the ice thickness distribution of all glaciers on earth. Nature Geoscience 12(3), 168. doi: 10.1038/s41561-019-0300-3.CrossRefGoogle Scholar
Favier, V, Wagnon, P and Ribstein, P (2004) Glaciers of the outer and inner tropics: a different behaviour but a common response to climatic forcing. Geophysical Research Letters 31(16), L16403. doi: 10.1029/2004GL020654.CrossRefGoogle Scholar
Fischer, M, Huss, M and Hoelzle, M (2015) Surface elevation and mass changes of all Swiss glaciers 1980–2010. The Cryosphere 9(2), 525540. doi: 10.5194/tc-9-525-2015.CrossRefGoogle Scholar
Fountain, AG, Krimmel, RM and Trabant, D and Geological Survey (U.S.) (1997) A Strategy for Monitoring Glaciers. Denver, CO: U.S. G.P.O.; Free on applications to the U.S. Geological Survey, Information Services, [Washington].CrossRefGoogle Scholar
Francou, B, Vuille, M, Favier, V and Cáceres, B (2004). New evidence for an ENSO impact on low-latitude glaciers: Antizana 15, Andes of Ecuador, 0°28′S. Journal of Geophysical Research: Atmospheres 109(D18). doi: 10.1029/2003JD004484.CrossRefGoogle Scholar
Francou, B, Vuille, M, Wagnon, P, Mendoza, J and Sicart, JE (2003) Tropical climate change recorded by a glacier in the central Andes during the last decades of the twentieth century: Chacaltaya, Bolivia, 16°S. Journal of Geophysical Research: Atmospheres 108(D5). doi: 10.1029/2002JD002959.CrossRefGoogle Scholar
Frans, C and 7 others (2015). Predicting glacio-hydrologic change in the headwaters of the Zongo River, Cordillera Real, Bolivia. Water Resources Research 51(11), 90299052. doi: 10.1002/2014WR016728.CrossRefGoogle Scholar
Gardelle, J, Berthier, E and Arnaud, Y (2012) Slight mass gain of Karakoram glaciers in the early twenty-first century. Nature Geoscience 5(5), 322325. doi: 10.1038/ngeo1450.CrossRefGoogle Scholar
Gilbert, A, Wagnon, P, Vincent, C, Ginot, P and Funk, M (2010) Atmospheric warming at a high-elevation tropical site revealed by englacial temperatures at Illimani, Bolivia (6340 m above sea level, 16°S, 67°W). Journal of Geophysical Research: Atmospheres 115(D10). doi: 10.1029/2009JD012961.CrossRefGoogle Scholar
Hanshaw, MN and Bookhagen, B (2014) Glacial areas, lake areas, and snow lines from 1975 to 2012: status of the Cordillera Vilcanota, including the Quelccaya Ice Cap, northern central Andes, Peru. The Cryosphere 8(2), 359376. doi: 10.5194/tc-8-359-2014.CrossRefGoogle Scholar
Haylock, M.R. and 23 others (2006) Trends in total and extreme South American rainfall in 1960–2000 and links with sea surface temperature. Journal of Climate 19(8), 14901512. doi: 10.1175/JCLI3695.1.CrossRefGoogle Scholar
Huss, M (2013) Density assumptions for converting geodetic glacier volume change to mass change. The Cryosphere 7(3), 877887. doi: 10.5194/tc-7-877-2013.CrossRefGoogle Scholar
Huss, M and Hock, R (2018). Global-scale hydrological response to future glacier mass loss. Nature Climate Change 8(2), 135. doi: 10.1038/s41558-017-0049-x.CrossRefGoogle Scholar
IPCC ed. (2014) Climate Change 2013 – The Physical Science Basis: Working Group I Contribution to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge: Cambridge University Press. Available at http://ebooks.cambridge.org/ref/id/CBO9781107415324.Google Scholar
Jordan, E (1991) Die Gletscher der Bolivianischen Anden. Steiner, Stuttgart.Google Scholar
Kääb, A, Berthier, E, Nuth, C, Gardelle, J and Arnaud, Y (2012) Contrasting patterns of early twenty-first-century glacier mass change in the Himalayas. Nature 488(7412), 495498. doi: 10.1038/nature11324.CrossRefGoogle ScholarPubMed
Kaser, G (2001) Glacier-climate interaction at low latitudes. Journal of Glaciology 47(157), 195204. doi: 10.3189/172756501781832296.CrossRefGoogle Scholar
Kaser, G and Osmaston, H (2002) Tropical Glaciers. Cambridge University Press, Cambridge.Google Scholar
Lejeune, Y and 7 others (2007) Melting of snow cover in a tropical mountain environment in Bolivia: processes and modeling. Journal of Hydrometeorology 8(4), 922937. doi: 10.1175/JHM590.1.CrossRefGoogle Scholar
Lejeune, Y (2009) Apports des modèles de neige CROCUS et de sol ISBAà l'etude du bilan glaciologique d'un glacier tropical et du bilan hydrologique de son bassin versant (PhD thesis), University of Grenoble, France.Google Scholar
Malz, P and 5 others (2018) Elevation and mass changes of the Southern Patagonia icefield derived from TanDEM-X and SRTM data. Remote Sensing 10(2), 188. doi: 10.3390/rs10020188.CrossRefGoogle Scholar
Mark, BG and Seltzer, GO (2005) Evaluation of recent glacier recession in the Cordillera Blanca, Peru (AD 1962–1999): spatial distribution of mass loss and climatic forcing. Quaternary Science Reviews 24(20–21), 22652280. doi: 10.1016/j.quascirev.2005.01.003.CrossRefGoogle Scholar
Maussion, F, Gurgiser, W, Großhauser, M, Kaser, G and Marzeion, B (2015) ENSO influence on surface energy and mass balance at Shallap Glacier, Cordillera Blanca, Peru. The Cryosphere 9(4), 16631683. doi: https://doi.org/10.5194/tc-9-1663-2015.CrossRefGoogle Scholar
McNabb, R, Nuth, C, Kääb, A and Girod, L (2019) Sensitivity of glacier volume change estimation to DEM void interpolation. The Cryosphere 13(3), 895910. doi: 10.5194/tc-13-895-2019.CrossRefGoogle Scholar
Morizawa, K, Asaoka, Y, Kazama, S and Gunawardhana, LN (2013) Temporal glacier area changes correlated with the El Niño/La Niña Southern Oscillation using satellite imagery. Hydrological Research Letters 7(2), 1822. doi: 10.3178/hrl.7.18.CrossRefGoogle Scholar
NASA JPL (2013) NASA Shuttle Radar Topography Mission Global 1 arc second Version 3. NASA JPL. doi: 10.5067/MEaSUREs/SRTM/SRTMGL1.003.Google Scholar
Nuth, C and Kääb, A (2011) Co-registration and bias corrections of satellite elevation data sets for quantifying glacier thickness change. The Cryosphere 5(1), 271290. doi: 10.5194/tc-5-271-2011.CrossRefGoogle Scholar
Paul, F and 19 others (2013) On the accuracy of glacier outlines derived from remote-sensing data. Annals of Glaciology 54(63), 171182. doi: 10.3189/2013AoG63A296.CrossRefGoogle Scholar
Rabatel, A and 7 others (2012) Can the snowline be used as an indicator of the equilibrium line and mass balance for glaciers in the outer tropics? Journal of Glaciology 58(212), 10271036. doi: 10.3189/2012JoG12J027.CrossRefGoogle Scholar
Rabatel, A and 27 others (2013) Current state of glaciers in the tropical Andes: a multi-century perspective on glacier evolution and climate change. The Cryosphere 7(1), 81102. doi: 10.5194/tc-7-81-2013.CrossRefGoogle Scholar
Rabatel, A, Machaca, A, Francou, B and Vincent, J (2006) Glacier recession on Cerro Charquini (16 S), Bolivia, since the maximum of the Little Ice Age (17th century). Journal of Glaciology 52(176), 110118. doi: 10.3189/172756506781828917.CrossRefGoogle Scholar
Ramirez, E and 8 others (2001) Small glaciers disappearing in the tropical Andes: a case-study in Bolivia: Glaciar Chacaltaya (16 o S). Journal of Glaciology 47(157), 187194.CrossRefGoogle Scholar
Réveillet, M, Rabatel, A, Gillet-Chaulet, F and Soruco, A (2015) Simulations of changes to Glaciar Zongo, Bolivia (16°S), over the 21st century using a 3-D full-Stokes model and CMIP5 climate projections. Annals of Glaciology 56(70), 8997. doi: 10.3189/2015AoG70A113.CrossRefGoogle Scholar
RGI Consortium (2017) Randolph Glacier Inventory – A Dataset of Global Glacier Outlines: Version 6.0: Technical Report, Global Land Ice Measurements from Space. NSIDC, Boulder. doi: 10.7265/N5-RGI-60.Google Scholar
Rolstad, C, Haug, T and Denby, B (2009) Spatially integrated geodetic glacier mass balance and its uncertainty based on geostatistical analysis: application to the western Svartisen ice cap, Norway. Journal of Glaciology 55(192), 666680. doi: 10.3189/002214309789470950.CrossRefGoogle Scholar
Sagredo, EA and Lowell, TV (2012) Climatology of Andean glaciers: a framework to understand glacier response to climate change. Global and Planetary Change 86–87, 101109.CrossRefGoogle Scholar
Schoolmeester, T (2018) The Andean Glacier and Water Atlas – The Impact of Glacier Retreat on Water Resources. UNESCO and GRID-Arendal, Paris, France and Arendal, Norway.Google Scholar
Seehaus, TC and 6 others (2016) Dynamic response of Sjögren inlet glaciers, Antarctic Peninsula, to ice shelf breakup derived from multi-mission remote sensing time series. Frontiers of Earth Science 4. doi: 10.3389/feart.2016.00066.Google Scholar
Seehaus, T and 5 others (2019) Changes of the tropical glaciers throughout Peru between 2000 and 2016 – mass balance and area fluctuations. The Cryosphere 13(10), 25372556. doi: 10.5194/tc-13-2537-2019.CrossRefGoogle Scholar
Sicart, JE, Hock, R, Ribstein, P, Litt, M and Ramirez, E (2011) Analysis of seasonal variations in mass balance and meltwater discharge of the tropical Zongo Glacier by application of a distributed energy balance model. Journal of Geophysical Research: Atmospheres 116(D13). doi: 10.1029/2010JD015105.CrossRefGoogle Scholar
Sicart, JE, Wagnon, P and Ribstein, P (2005) Atmospheric controls of the heat balance of Zongo Glacier (16°S, Bolivia). Journal of Geophysical Research: Atmospheres 110(D12). doi: 10.1029/2004JD005732.CrossRefGoogle Scholar
Silverio, W and Jaquet, J-M (2017) Evaluating glacier fluctuations in Cordillera Blanca (Peru). Archives Des Sciences/Editees Par La Societe De Physique Et D'histoire Naturelle De Geneve 18, 145162.Google Scholar
Soruco, A and 9 others (2009 b) Mass balance of Glaciar Zongo, Bolivia, between 1956 and 2006, using glaciological, hydrological and geodetic methods. Annals of Glaciology 50(50), 18. doi: 10.3189/172756409787769799.CrossRefGoogle Scholar
Soruco, A and 6 others (2015) Contribution of glacier runoff to water resources of La Paz city, Bolivia (16°S). Annals of Glaciology 56(70), 147154. doi: 10.3189/2015AoG70A001.CrossRefGoogle Scholar
Soruco, A, Vincent, C, Francou, B and Gonzalez, JF (2009 a) Glacier decline between 1963 and 2006 in the Cordillera real, Bolivia. Geophysical Research Letters 36(3), L03502. doi: 10.1029/2008GL036238.CrossRefGoogle Scholar
Toutin, T (2002) Three-dimensional topographic mapping with ASTER stereo data in rugged topography. IEEE Transactions on Geoscience and Remote Sensing 40(10), 22412247. doi: 10.1109/TGRS.2002.802878.CrossRefGoogle Scholar
Troll, C (1941) Studien zur vergleichenden Geographie der Hochgebirge der Erde. Verlag nicht ermittelbar, Bonn.Google Scholar
Veettil, BK, Wang, S, Simões, JC and Pereira, SFR (2018) Glacier monitoring in the eastern mountain ranges of Bolivia from 1975 to 2016 using Landsat and Sentinel-2 data. Environmental Earth Sciences 77(12), 452. doi: 10.1007/s12665-018-7640-y.CrossRefGoogle Scholar
Vijay, S and Braun, M (2016) Elevation change rates of glaciers in the Lahaul-Spiti (Western Himalaya, India) during 2000–2012 and 2012–2013. Remote Sensing 8(12), 1038. doi: 10.3390/rs8121038.CrossRefGoogle Scholar
Vijay, S and Braun, M (2018) Early 21st century spatially detailed elevation changes of Jammu and Kashmir glaciers (Karakoram–Himalaya). Global and Planetary Change 165, 137146. doi: 10.1016/j.gloplacha.2018.03.014.CrossRefGoogle Scholar
Vuille, M and 6 others (2008) Climate change and tropical Andean glaciers: past, present and future. Earth-Science Reviews 89(3–4), 7996. doi: 10.1016/j.earscirev.2008.04.002.CrossRefGoogle Scholar
Vuille, M and 12 others (2018) Rapid decline of snow and ice in the tropical Andes – impacts, uncertainties and challenges ahead. Earth-Science Reviews 176(Suppl. C), 195213. doi: 10.1016/j.earscirev.2017.09.019.CrossRefGoogle Scholar
Vuille, M and Bradley, RS (2000) Mean annual temperature trends and their vertical structure in the tropical Andes. Geophysical Research Letters 27(23), 38853888. doi: 10.1029/2000GL011871.CrossRefGoogle Scholar
Vuille, M, Bradley, RS, Werner, M and Keimig, F (2003) 20th century climate change in the tropical Andes: observations and model results. Climatic Change 59(1–2), 7599. doi: 10.1023/A:1024406427519.CrossRefGoogle Scholar
Wagnon, P, Ribstein, P, Francou, B and Pouyaud, B (1999) Annual cycle of energy balance of Zongo Glacier, Cordillera Real, Bolivia. Journal of Geophysical Research: Atmospheres 104(D4), 39073923. doi: 10.1029/1998JD200011.CrossRefGoogle Scholar
Wagnon, P, Ribstein, P, Francou, B and Sicart, JE (2001) Anomalous heat and mass budget of Glacier Zongo, Bolivia, during the 1997/98 El Niño year. Journal of Glaciology 47(156), 2128. doi: 10.3189/172756501781832593.CrossRefGoogle Scholar
WGMS. Fluctuations of Glaciers 1990–1995. WGMS.Google Scholar
Zemp, M and 14 others (2019) Global glacier mass changes and their contributions to sea-level rise from 1961 to 2016. Nature 1, 382386. doi: 10.1038/s41586-019-1071-0.CrossRefGoogle Scholar
Zink, M, Bartusch, M and Miller, D (2011) TanDEM-X Mission Status. IEEE International Geoscience and Remote Sensing Symposium (IGARSS), pp. 14. Available at http://elib.dlr.de/70355/.CrossRefGoogle Scholar
Figure 0

Fig. 1. (a) Map of Bolivia. Dashed square indicates the study area, blue polygons show the glacier extent based on the RGI 6.0. Background: Natural Earth (b) Topographic map of analyzed mountain ranges and TanDEM-X (TDX) and Pléiades data coverage. Background: SRTM © NASA.

Figure 1

Fig. 2. Area on (light blue) and off glacier (red) and off glacier dh/dt distribution (blue dots) based on slope. Error bars indicate NMAD of dh/dt values for each slope bin. Dotted line shows the slope threshold used in the analysis. Note: the glacier area is scaled by a factor of 10.

Figure 2

Fig. 3. Glacier area changes. Ice divides (black polygons) are from the glacier outlines in 2000. Dashed polygons indicate subsets illustrated in Figure S3 (Supplement) Background: SRTM hillshade © NASA 2000.

Figure 3

Fig. 4. Relative area changes (dS) (2000–2013) of individual glaciers (dot color) plotted against glacier size (dot size), median elevation (distance from center) and mean aspect (orientation). Red circle: equilibrium line altitude (ELA) from Rabatel and others (2012).

Figure 4

Table 1. Observed glacier changes for different time intervals

Figure 5

Fig. 5. Glacier surface elevation changes (unfiltered; left panel: 2000–2013; right panel: 2013–2016). Dashed polygons indicate subsets illustrated in Figures S6 and S7 (Supplement). Background SRTM hillshade © NASA 2000.

Figure 6

Fig. 6. Hypsometric distribution of glacier area (light blue), glacier area with dh/dt measurements (red) and mean, filtered dh/dt values (blue dots) of each hypsometric bin for the observation period 2000–2013. Error bars indicate NMAD of dh/dt for each hypsometric bin.

Figure 7

Fig. 7. Specific mass balance (spMB) (2000–2013) of individual glaciers (dot color) plotted against glacier size (dot size), median elevation (distance from center) and mean aspect (orientation). Red circle: equilibrium line altitude (ELA) from Rabatel and others (2012). Note: only glaciers with >40% elevation change data coverage, which is spread over >2/3 of the hypsometric distribution are included.

Figure 8

Fig. 8. Temporal evolution of glacier area and specific mass balance (spMB) (a and b) as well as Oceanic Nino Index (ONI, c) for the studied period.

Supplementary material: PDF

Thorsten et al. supplementary material

Thorsten et al. supplementary material Figures S1-S20 and Table S1

Download Thorsten et al. supplementary material(PDF)
PDF 4 MB
You have Access Open access
9
Cited by

Save article to Kindle

To save this article to your Kindle, first ensure coreplatform@cambridge.org is added to your Approved Personal Document E-mail List under your Personal Document Settings on the Manage Your Content and Devices page of your Amazon account. Then enter the ‘name’ part of your Kindle email address below. Find out more about saving to your Kindle.

Note you can select to save to either the @free.kindle.com or @kindle.com variations. ‘@free.kindle.com’ emails are free but can only be saved to your device when it is connected to wi-fi. ‘@kindle.com’ emails can be delivered even when you are not connected to wi-fi, but note that service fees apply.

Find out more about the Kindle Personal Document Service.

Mass balance and area changes of glaciers in the Cordillera Real and Tres Cruces, Bolivia, between 2000 and 2016
Available formats
×

Save article to Dropbox

To save this article to your Dropbox account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you used this feature, you will be asked to authorise Cambridge Core to connect with your Dropbox account. Find out more about saving content to Dropbox.

Mass balance and area changes of glaciers in the Cordillera Real and Tres Cruces, Bolivia, between 2000 and 2016
Available formats
×

Save article to Google Drive

To save this article to your Google Drive account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you used this feature, you will be asked to authorise Cambridge Core to connect with your Google Drive account. Find out more about saving content to Google Drive.

Mass balance and area changes of glaciers in the Cordillera Real and Tres Cruces, Bolivia, between 2000 and 2016
Available formats
×
×

Reply to: Submit a response

Please enter your response.

Your details

Please enter a valid email address.

Conflicting interests

Do you have any conflicting interests? *