Mass Balance of Four Mongolian Glaciers: In-situ Measurements, Long-Term Reconstruction and Sensitivity Analysis

This study describes the in situ observed and modeled mass balances of four representative glaciers in the Mongolian Altai Mountains. Mass-balance stakes and automatic weather stations (AWS) were installed along each glacier to obtain the in situ mass balances and meteorological variables. We calibrate the ERA5 meteorological variables using the observed values, then estimate the precipitation parameters to obtain the observed point mass balances using an energy- and mass-balance model. We evaluate the mass-balance profile and glacier-wide mass balance via a mass-balance model using the calibrated ERA5 data and estimated precipitation parameters to fill in the spatial and temporal gaps in the stake networks and AWS measurements. We demonstrate that almost all of the observed mass balances have been in a negative state since 2003. We also reconstruct the long-term mass balances for the 1980–2018 period, which range from −760 to −160 mm water equivalent under the different climatic conditions (annual precipitation varies from 190 to 860 mm). We further evaluate the mass-balance sensitivities to temperature and precipitation changes. We confirm that both sensitivities correlate significantly with the annual precipitation; increased precipitation yields more negative sensitivity to temperature changes and less positive sensitivity to precipitation changes.


INTRODUCTION
Mongolia is a land-locked arid country in Northeast Asia (inset map of Figure 1A) whose isolated location has resulted in limited moisture flux from the ocean (Sato et al., 2007;Antokhina et al., 2019). However, the regional climate, especially within the glacierized zone, remains equivocal. Furthermore, the impact of climate change on the high-mountain landscape, including glaciers, is poorly understood (Ganyushkin et al., 2017). Mongolian glaciers form within two major mountain ranges, the Altai and Khangai mountains, both of which extend from northwest to southeast. Mongolia generally receives limited precipitation because of its arid to semi-arid climate, with 300-400 mm of mean annual precipitation in the Khangai Mountains and 250-300 mm in the Altai Mountains (Dagvadorj et al., 2009). Approximately 10% of the freshwater volume in Mongolia is estimated to be stored as glacier ice (Dashdeleg et al., 1983;Baast, 1998;Myagmarjav and Davaa, 1999), with many of the river basins already experiencing water stress (Walther et al., 2017). Pan et al. (2019) estimated that about 10% of river runoff in the Ulgii catchment was attributed to glacier meltwater. The spatial glacier distribution is sporadic and decreases from northwest to southeast ( Figure 1A). They are categorized as summeraccumulation-type glaciers (Supplementary Figure S1), which indicates that the glacier mass balances are particularly sensitive to temperature changes (Fujita and Ageta, 2000;Fujita, 2008). Most of the glaciers have been retreating recently, which raises the concern that glacier meltwater will not continue to make a sustainable contribution to the regional water supply in the future (Davaa et al., 2006;Kadota and Davaa, 2007;Kadota et al., 2011).
The first evidence of changes on Mongolian glaciers was reported by Potanin (1883) after his expeditions to the Mongolian Altai (1877-1879) and a subsequent Polish geological expedition in 1966 (Rutkowski and Slowanski, 1970). In 1910, a Royal Geographical Society expedition conducted an extensive survey of the mountains, produced a detailed topographic map, and documented the glacier extents based on photographs from the Turgen Mountains (Carruthers, 1912). A United States-Mongolian expedition retraced portions of the 1910 expedition in summer 2010, and analyzed repeat field observations and photographs that were acquired in 1910 and 2010 (Kamp et al., 2013). The lower portions of the glaciers typically underwent significant retreat, whereas the ice volumes in the glacier heads appeared to be stable. For example, West Turgen Glacier retreated by ∼600 m and downwasted by ∼70 m during the 1910-2010 period (Kamp et al., 2013). Landsat-based regional glacier inventories have been created from the 1990, 2000, and 2010 imagery (Kamp and Pan, 2015), and glacier recession has been analyzed in detail (Pan et al., 2018). However, repeated in situ observations were limited until the 2000s.
The Tsambagarav and Potanin glaciers have been monitored continuously since 2003 as part of a collaboration between the Institute of Meteorology and Hydrology, Mongolia, and the Japan Agency for Marine-Earth Science and Technology (JAMSTEC; Kadota et al., 2011). The mass balance of Potanin Glacier has been determined using glaciological and meteorological observations (Konya et al., 2010), with similar mass-balance fluctuations to those observed along Maliy Aktru Glacier in the Russian Altai Mountains (Konya et al., 2013). Nakazawa et al. (2015) analyzed the pollen concentrations and oxygen stable isotopes in a snow pit  Sakai, 2019). Black circles and blue rhombus in (A) denote the meteorological stations, whose climate data were used in Supplementary Figure S1, and Maliy Aktru Glacier, whose observed mass balance was compared with the estimation of this study ( Figure 8C and Supplementary Table S6), respectively. Red polygons denote the glacier outlines, which have been delineated from the Sentinel-2A images (Table 1). Black dots and blue triangles denote the glacier mass-balance stakes and automatic weather stations (AWSs), respectively. Black cross on Potanin Glacier denotes a pit site at which accumulation was analyzed by pollen analysis (Nakazawa et al., 2015). on Potanin Glacier and estimated the snow-accumulation and melt-intensity histories during the 2007-2011 period. The Mongolian Government's Information and Research Institute of Meteorology, Hydrology, and Environment has continued these glacier-monitoring studies and expanded their efforts to acquire glaciological observations on other Mongolian glaciers since JAMSTEC reduced its involvement in 2013. Although a growing wealth of in situ glaciological observations have been acquired and satellite-based aerial changes revealed (Pan et al., 2018), relatively little is known about the fluctuations in glacier mass across the continental climates of Mongolia. This study aims to summarize the glaciological observations since 2003 by estimating the precipitation regimes of four Mongolian glaciers using an energy-and mass-balance model, reconstruct their longterm mass balances , and reveal their mass-balance sensitivities.

Glaciers
Four Mongolian glaciers were selected for this mass-balance study (Figure 1). Potanin and Turgen are valley-type glaciers, whereas Tsambagarav and Sutai are ice-cap-like glaciers; the location, elevation range, area and first year of observations for each glacier are summarized in Table 1. We determined the glacier hypsometry at 50-m elevation intervals using the manually delineated glacier outlines from 2017 satellite images ( Figure 1) and ASTER-GDEM2 (Tachikawa et al., 2011). Four automatic weather stations (AWSs) were installed near the glaciers in 2005 (Tsambagarav), 2012 and 2014 (Potanin), and 2013 (Turgen) to monitor local meteorological parameters ( Table 2).
Potanin Glacier is located in the Tavan Bogd region, which is along the border between Mongolia, Russia, and China ( Figure 1A). The Tavan Bogd region is designated as a national park, thereby protecting it from residential development and making it a popular tourist destination. Potanin Glacier flows from the summit of Mt. Khuiten (4,374 m above sea level (a.s.l.)), which is the highest peak in the Altai Mountains, toward the east, extending 10.4 km along a valley and terminating at 2,907 m a.s.l. ( Figure 1B). The glacier area is estimated to be 24.7 km 2 .
Sutai Glacier is located on Sutai Mountain (4,217 m a.s.l.), the highest peak in the Gobi-Altai Mountains. The mountain range stretches 60 km along the border between the Khovd and Govi-Altai provinces ( Figure 1A). The glacier flows from the summit of Sutai Mountain to 3,459 m a.s.l, and possesses an estimated glacier area of 6.2 km 2 ( Figure 1C). Although the glacier is identified as multiple (14) glaciers in the GAMDAM glacier inventory (Sakai, 2019), we deal with this as a single glacier.
Tsambagarav Glacier is located in the Tsambagarav Mountains, which are part of the Altai Mountains and lie along the border between the Khovd and Bayan-Ulgii provinces in western Mongolia ( Figure 1A). The glacier flows 3.6 km, from 3,814 to 3,226 m a.s.l, and possesses an estimated glacier area of 6.3 km 2 ( Figure 1D). Although the glacier is identified as multiple (4) glaciers in the GAMDAM glacier inventory (Sakai, 2019), we deal with this as a single glacier.  Figure 1 Aug. 10, 2017Aug. 11, 2017Sep. 18, 2017Sep. 18, 2017 zPR denotes the reference elevation for calculating precipitation gradient. Turgen Glacier is located in the Turgen Mountains, which are in northwestern Mongolia, ∼80 km south of the Russian border ( Figure 1A). The mountains are heavily glacierized on their northern slopes and feed several rivers that flow northeast toward the regional capital of Ulaangom and into the Lake Uvs Nuur Basin. The glacier flows from south (3,905 m a.s.l.) to north (3,020 m a.s.l.), extending for 3.6 km, and possesses an estimated glacier area of 3.3 km 2 .

Mass-Balance Observations
We installed mass-balance stakes at 50-100-m elevation intervals on the Potanin (14 stakes), Tsambagarav (5 stakes), Turgen (7 stakes), and Sutai (4 stakes) glaciers, and acquired almost annual measurements since 2003, 2010, 2013, and 2015, respectively ( Figure 1). We measured the stake height from the ice surface and reinstalled the stakes when necessary. The stake positions were measured intermittently using a differential global positioning system to obtain the flow velocity (not shown in this study). The stakes were sometimes moved up-glacier to compensate for the glacier flow. The stake height were converted into water equivalent (w.e.) with an ice density of 870 kg m −3 . In the most measurements, snow thickness above the ice surface was zero or a few centimeter so that we ignored snow in the mass balance estimation.

Mass-Balance Model
The stakes that were installed for the mass-balance measurements did not cover the entire glacier. The stakes on Potanin and Turgen glaciers primarily covered the lower portions of the glaciers, whereas those on Sutai and Tsambagarav glaciers primarily covered the upper portions ( Figure 1). We calibrated the precipitation parameters using an energy-and mass-balance model, which has been used in previous mass-balance studies of various Asian glaciers (e.g., Fujita and Ageta, 2000;Sakai et al., 2009;Sunako et al., 2019), to evaluate the glacier-wide mass balances from these limited observations. We also reconstructed their long-term mass balances  to understand how these glaciers have responded to recent climate change.
The model used in this study calculated the daily energy balance as: where Q M is the excess energy (heat) for melt (W m −2 ); α is the surface albedo (dimensionless); R S is the downward solar radiation (W m −2 ); R L is the downward long-wave radiation (W m −2 ); ε is the emissivity of the snow/ice surface, which is assumed to be 1.0 (dimensionless); σ is the Stefan-Boltzmann constant (5.67 × 10 −8 W m −2 K −4 ); T S is the surface temperature (°C); H S is the sensible heat (W m −2 ); H L is the latent heat (W m −2 ); and Q G is the heat that is conducted into the snow and ice (W m −2 ). All of the down-going components have a positive sign. The excess energy (Q M ) is used to melt the snow/ice at the surface once the surface temperature reaches 0°C, and its daily amount is obtained as: where M S is the daily amount of snow/ice melt (kg m −2 d −1 or mm w. e. d −1 ); t d is the length of a day in seconds (86,400 s); and l M is the latent heat for ice melt (3.33 × 10 5 J kg −1 ). The turbulent heat fluxes (H S and H L ) are estimated via a bulk aerodynamic method (Kondo, 1994;Fujita and Sakai, 2014). The surface albedo (α) is estimated via a scheme that models the exponential decay of the snow albedo after a fresh snowfall (Kondo and Xu, 1997;Fujita and Sakai, 2014). The ice albedo, which is a constant parameter in the model, was assumed to be 0.25 based on an observation in 2019 (Supplementary Figure  S2). The downward long-wave radiation (R L ) is estimated from the air temperature, relative humidity, and solar radiation, and is used as a proxy for cloudiness based on the empirical equations proposed by Kondo (1994). We compared the estimated downward long-wave radiation with the observations from the Potanin High AWS ( Table 2) and established a calibration equation for the studied glaciers. Precipitation is distinguished as either solid snow or liquid rain, as follows: where r sr is the snowfall probability, which changes linearly from 1 to 0 (dimensionless); T d , T sn , and T rn are the daily air temperature, and threshold temperatures for snow and rain (°C), respectively; and P S and P d are the daily amounts of snowfall and precipitation (mm w. e. d −1 ), respectively. The threshold temperatures for snow and rain are assumed to be 0 and 8°C, respectively, based on in situ observations by Konya et al. (2013). The specific mass balance at a given site or elevation (b z , m w. e.) is then obtained as: where E V is the daily amount of evaporation (mm w. e. d −1 ), which is converted from the latent heat (H L ); R F is the daily amount of refreezing water (mm w. e. d −1 ); and ρ W is the water density (kg m −3 ) for the unit conversion (mm to m). The daily mass balance is summed for a given period (observation period or 1 year). The glacier-wide mass balance (B, m w. e.) is calculated using the hypsometry as follows: where A z and b z are the glacier area and mass balance for a given 50-m elevation band, respectively. Detailed descriptions of the model are available in Fujita and Ageta (2000) and Fujita and Sakai (2014). The model has been applied to estimate the mass Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 balance and runoff of glaciers in an arid region of China, which possesses a similar environment to that of the glaciers in this study (Sakai et al., 2009, Sakai et al., 2010.

Reanalysis Data
We used the ERA5 reanalysis dataset (0.25°× 0.25°), which is produced by the European Centre for Medium-Range Weather Forecasts (Hersbach et al., 2020), as the input meteorological variables for the mass-balance model to estimate the glacierwide mass balance, reconstruct the long-term mass balance, and evaluate the glacier mass-balance sensitivity for each of the studied glaciers. We tested three methods for the air temperature extraction from the reanalysis data. We extracted the daily air temperature at each AWS site (T PL ,°C ) from the daily pressure level temperatures (T p and T p+1 ,°C ) at the closest geopotential heights (z p and z p+1 , m a.s.l.) containing the site elevation (z site , m a.s.l.) via the first method as: where p is a pressure level that satisfies z p ≤ z site < z p+1 , and L PL is the temperature lapse rate (°C m −1 ). We modified the daily 2-m height temperature (T 2m ,°C) using the orographic height (taken as the surface elevation in the ERA5 reanalysis; z ERA5 , m a.s.l.), site elevation, and lapse rate between the two elevations via the second method as: where p u and p l are the pressure levels that satisfy z p u ≤ z site < z p u +1 and z p l ≤ z ERA5 < z p l +1 , with the general condition that z ERA5 < z site , respectively; and L C is the temperature lapse rate (°C m −1 ). We also tested a third method for the temperature extraction to avoid the use of hourly ERA5 pressure level data, since these data require a huge amount of storage. The monthly temperature at the site elevation (T PL ,°C) was estimated from the monthly pressure level temperatures with the same manner of Eq. 6. The daily air temperature was then estimated using the daily 2-m temperature (T 2m ,°C) and its monthly mean (T 2m ,°C) as: We used these estimated air temperatures to compare the daily parameters (solar radiation, wind speed, and relative humidity) of the ERA5 products with the observational data at the Potanin, Tsambagarav, and Turgen AWSs. We adopted the pressure-levelbased temperature (T PL ), then obtained a linear regression equation for each parameter (see Section 3.1 for the full details). For Sutai Glacier, we assumed that the terminus elevation represented the elevation for the temperature extraction (3,818 m a.s.l.), as no AWS data were available due to battery issues.

Precipitation Calibration
A tipping bucket rain gauge, which is commonly installed at AWSs, cannot capture the precipitation amount precisely because snow is the dominant precipitation type, even during the summer melt season. Therefore, we calibrated the precipitation parameter in this study by applying a multiplier to the ERA5 daily precipitation (r p , dimensionless). We determined the precipitation parameter that yielded the same value as the observed mass balance for each stake measurement period (measurements at a given stake for a given period), and then obtained the average and standard deviation for each measurement period. The daily precipitation (P z , mm w. e.) at a given elevation (z, m a.s.l.) was then calculated as: where P e is the ERA5 daily precipitation (mm w. e.); L p is the precipitation gradient (% km −1 ), which was assumed to be 45% km −1 based on the estimated precipitation parameter for Potanin Glacier (Section 3.2); and z PR is the reference elevation for estimating the precipitation along the given elevation, which was obtained from the mean elevation of the mass-balance stakes across each glacier. The denominator value (100,000) is for the unit conversion (% km −1 to m −1 ). We also determined a single precipitation parameter for each glacier to reconstruct the long-term mass balance, and then verified the reliability of each parameter by calculating the root mean square error (RMSE) and mean error (ME) between the estimated and measured point mass balances.

Mass-Balance Reconstruction, Climate Regime, and Sensitivity
We reconstructed the long-term annual mass balances for the 1980-2018 period using the determined precipitation parameter for each glacier to better understand how air temperature and precipitation variations can affect glacier mass balance. We also determined the climatic regime for the present-day glacier geometry to demonstrate how the current climatic conditions affect the glacier mass balance by varying the air temperature and precipitation amount (Fujita et al., 2011;Sunako et al., 2019). The daily air temperature (T d ,°C) and precipitation (P d , mm) for a given year are defined as: where T and P denote air temperature (°C) and precipitation (mm), respectively; the d and e subscripts denote the daily variables for the model input and ERA5, respectively; and the E and C subscripts denote the mean summer T and annual sum of P from the ERA5 data and controlled variable, respectively. We varied the controlled summer mean temperature (T C ) between −4 and +4°C (1°C interval), and the annual precipitation (P C ) between 40 and 1,600 mm (20-mm interval), to obtain the range of glacier-wide mass balances (B) for a given year of ERA5 data.
Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 We finally averaged the glacier-wide mass balances for the 39years ERA5 patterns. We performed a sensitivity analysis of the glacier mass balance to understand how the four glaciers would respond to temperature (±0.5°C) and precipitation (±5%) changes. We also evaluated the uncertainties in the estimated glacier-wide mass balances using the precipitation parameter (±1σ) and precipitation gradient (±45% km −1 ) uncertainties.

Meteorological Data
We first compared the daily ERA5 air temperatures that were calibrated via the three above-mentioned methods (Section 2.4) with the in situ AWS temperatures at the Potanin, Tsambagarav, and Turgen glaciers (Supplementary Figure S3). The ERA5 temperatures that were calibrated via methods 1 (T PL ) and 3 (T ML ), which are based on the pressure level temperature, are highly consistent with the observed temperatures (T obs ) in terms of the slope (0.99-1.04) and R 2 (0.90-0.98) of the regression lines, and RMSE (1.3-3.2°C), whereas the method 2 calibration introduced a significant cold bias during the cold season (center panels of Supplementary Figure S3). Method 2 (T C2m ) calibrates the 2-m temperature in the ERA5 dataset by considering the temperature lapse rate between the site elevation and ERA5 orography. A strong inversion layer usually develops across Mongolia during winter (Dagvadorj et al., 2009). The 2-m temperature in the ERA5 dataset could represent such a cold bias near the winter surface, whereas the glaciers could be situated above this inversion layer. We therefore adopted the calibrated ERA5 air temperature via method 1 (T PL ) for the following calculations in this study.
We also compared the solar radiation, wind speed, and relative humidity in the reanalysis dataset with the in situ observations to further examine the performance of the ERA5 dataset for each glacier (Supplementary Figures S4-S6). Strong correlations are observed for the daily solar radiation (R 2 0.77-0.88). However, only very weak to moderate correlations are observed for the wind speed and relative humidity (R 2 0.10-0.55 for wind speed, and R 2 0.03-0.59 for relative humidity), even though the p-values indicate a degree of significance owing to the large sample numbers. The strong air temperature and solar radiation correlations are due to their obvious seasonality, such that they do not guarantee the reliability of the reanalysis data; however, the small RMSE values support the representativeness of the reanalysis variables. Therefore, we adopted the linear regression equations for air temperature and solar radiation to calibrate the ERA5 data for the massbalance calculations. We used the uncalibrated ERA5 data when a regression was unavailable (air temperature and solar radiation at Sutai Glacier and solar radiation at Turgen Glacier). The ERA5 data clearly underestimated the wind speed ( Supplementary  Figures S4C-S6C). Therefore, we simply assumed the wind speed as u 3.0 u ERA (where u and u ERA are the model input and ERA5 wind speeds, respectively) based on the relationships shown in Supplementary Figures S4C-S6C. We calibrated the ERA5 humidity data by subtracting 0.11 based on the obtained relationships for the Potanin and Turgen glaciers (Supplementary Figure S4D and Supplementary Figure S6D) We checked the performances of both the ERA5 air temperature lapse rate and modeled long-wave radiation against the Potanin Glacier observations (Supplementary Figure S7). The temporal changes in the air temperature lapse rate exhibit similar seasonal patterns, with less negative changes during winter and more negative changes during summer. (Supplementary Figure S7A). The less negative lapse rate during winter is likely due to strong inversion conditions, as indicated by the calibrated ERA5 temperatures. Conversely, although the temperature lapse rate along the glacier surface tends to be small (less negative) due to the katabatic-wind-driven cooling effect (Shea and Moore, 2010;Radić et al., 2014), the temperature lapse rate at Potanin Glacier is rather large (more negative) during summer, which is likely due to adiabatic dry conditions. The long-wave radiation, which was estimated from the air temperature, relative humidity, and solar radiation as a proxy for cloud cover (Kondo, 1994;Fujita and Ageta, 2000), also exhibited similar temporal fluctuations to the in situ observations (Supplementary Figure S7C). Scatter and temporal plots of the temperature lapse rate and long-wave radiation results across Potanin Glacier suggest that the linearly calibrated variables (LAP calib and LRD calib ) are good representatives of the observed values. We therefore adopted this linear regression approach to obtain the air temperature lapse rate and longwave radiation values for the other glaciers in this study.

Determining the Precipitation Parameters
We used the calibrated ERA5 data to determine the precipitation parameter that minimized the misfit between the modeled and observed mass balances via iterative calculations at each stake measurement ( Figure 2). Histograms of the precipitation parameter indicate that there are reasonable distributions for the Potanin,Tsambagarav,and Turgen glaciers (Figures 2A,C,D), whereas that for Sutai Glacier is rather flat ( Figure 2B); this high degree of variability is likely due to the small number of samples compared with those for the other three glaciers (vertical width of the box plot corresponds to the sample number).
We estimated the precipitation gradient to calculate the mass balance at the unmeasured elevations (Supplementary Figure  S8) using pollen-based accumulation rates (Nakazawa et al., 2015). The slope of regression line suggests that the precipitation increased by 44% per kilometer of elevation; we therefore assumed a precipitation gradient (L p ) of 45% km −1 for all four glaciers. We then calculated the mass-balance profile for each measurement period (Figures 3-6). The linear regressions were obtained for band-mean specific mass balances, which were averaged at a 50-m elevation band. The mass balance that was obtained via a linear regression of the stake measurements along Potanin Glacier tended to be more positive across the highelevation zone, with the exception of the first measurements during the 2003-2004 period (Figure 3). The precipitation Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 parameters were estimated to be less than one in many cases. The linear regression of the stake measurements is less than the simulated profile along the lower elevations of Sutai Glacier, with the exception of the first measurements during the 2015-2016 period (Figure 4). The estimated precipitation parameters vary significantly among the four periods, resulting in the flat distribution in the histogram ( Figure 2B). There were two periods (2010-2011 and 2016-2017) that exhibited unusual linear regressions for Tsambagarav Glacier ( Figure 5), with three periods exhibiting similar profiles between the linear regression and simulation. The modeled profile for Turgen Glacier was much more positive than that obtained via the linear regression approach (Figure 6). These results and the obtained glacier-wide mass balances are summarized in Supplementary Tables S2-S5. The estimated precipitation parameters for each stake measurement were averaged to determine a single long-term parameter value for each glacier (Figure 2 and Table 3). The determined precipitation parameter (r p ) for each glacier yielded the best estimate of the stake observations, whereby the RMSEs were minimized and the MEs were slightly positive. The ERA5 precipitation data (P ERA5 ) were in the 300-840 mm range for the four glaciers, with the maximum occurring across Sutai Glacier (Table 3). However, the calibrated annual precipitation value (P cal ) for Turgen Glacier had increased slightly to 860 mm, with this increase likely due to the precipitation gradient. The precipitation values for the other glaciers were significantly reduced, ranging between 190 and 640 mm ( Table 3).

Long-Term Reconstruction
The determined precipitation parameters have been used to calculate the long-term mass-balance profiles averaged for the 1980-2018 period. Figure 7 shows the hypsometries (50-m elevation bands) and mass-balance profiles for the four studied glaciers. The Sutai and Tsambagarav glaciers, which are ice-caplike, have top-heavy hypsometries, and their mass-balance profiles show that their equilibrium line altitudes (ELAs) are situated along the uppermost portions of the glaciers, suggesting glacier-wide negative mass balances. Conversely, the ELAs of the Potanin and Turgen glaciers, which are valley-type, are situated within the current glacier geometry, and at lower elevations than those of the two ice-cap-like glaciers. The shaded regions along the thick solid lines denote the uncertainty due to the estimated precipitation parameters, and the dotted and dashed lines denote the mass-balance profiles with zero (dotted lines) and doubled (90% km −1 ) precipitation gradients.
The time series of the summer (June, July, and August) mean temperature anomaly (dSMT), calibrated annual precipitation, mass balance, and ELA for the four studied glaciers is shown in Figure 8. The summer mean temperatures (SMTs), which should strongly affect the ablation, are obtained by weighting the temperatures with the hypsometry ( Table 3). The SMT anomalies fluctuate in a similar manner among the glaciers and show a significant warming trend from the mid-1990s ( Figure 8A). However, there is no significant trend in precipitation among the four glaciers, with contrasting FIGURE 2 | Histogram and box plot (upper panels), and root mean square error (RMSE) and mean error (ME) (lower panels) of the precipitation parameter for Potanin (A, E), Sutai (B, F), Tsambagarav (C, G) and Turgen (D, H) glaciers. The horizontal box width, line, circle and whiskers denote the interquartile, median, mean and standard deviations, respectively. The vertical width of the box indicates the number of mass balance measurements used to determine the precipitation parameter (n). The vertical line and grey shaded region in the lower panels denote the mean and standard deviations of the precipitation parameter, respectively.
Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 fluctuations observed. Very arid conditions and similar fluctuations are observed for the Sutai and Tsambagarav glaciers, two ice-cap-like glaciers that are situated in the southeastern region, a moderate amount of precipitation is observed for Potanin Glacier, and a large amount of precipitation is observed for Turgen Glacier ( Figure 8B). The glacier-wide mass balances of the four glaciers are summarized as follows: Sutai Glacier has been in a negative balance for the entire analysis period (1980-2018, yellow), whereas the other three glaciers have been in a negative mass balance since the mid-1990s ( Figure 8C).  Figure 8D).
In Figure 8C, the in-situ observational mass balances of Maliy Aktru Glacier in Russian Altai (location shown in Figure 1A) and of Urumqi No. 1 Glacier in Tienshan (location not shown) are also depicted (WGMS, 2021). Interannual variabilities and longterm trends of the estimated and observed mass balances seem . The blue dots denote the observed mass balance, and the blue line denotes a linear regression obtained for band-mean specific mass balances, which were averaged at a 50-m elevation band, respectively. The black line and grey shaded region denote the modelled mass-balance profile and its variability due to the precipitation parameter (r p in each panel), respectively. The thin red line denotes the modelled mass-balance profile without the precipitation calibration. Orange dots denote the mass balance determined via a pollen analysis in a pit (Nakazawa et al., 2015). The correlation coefficients among the glacier-wide mass balance (B), area-weighted SMT, P cal , and ELA are summarized in Table 3. There are significant correlations among each of the variables (p < 0.001), with the SMTs showing stronger correlations than P cal . The weak to moderate negative correlations between SMT and P cal suggest that warm conditions tend to coincide with more arid conditions (Table 3; Figures 8A,B).

Climatic Regime and Sensitivity
The climatic regimes for the present-day geometries of the four glaciers as a function of the area-weighted SMT and annual  Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 9 precipitation are shown in Figure 9. The calibrated SMT and annual precipitation values yield a negative mass balance for the Sutai and Tsambagarav glaciers ( Figures 9B,C), whereas there is a large degree of scatter among the variables, which yield a positive mass balance for part of the Potanin and Turgen glaciers ( Figures 9A,D).
We further tested the mass-balance sensitivity by varying the air temperature (±0.5°C) and precipitation (±5%) during the simulations ( Table 3). The sensitivity to temperature changes is summarized by determining the resultant negative mass balance due to 1°C warming, which ranges between −500 and −250 mm w. e. (dB/dT). However, the sensitivity to precipitation changes (dB dP10% ) cannot be compared in such a simplistic manner because the degree of the precipitation change differs among both the glaciers and the years in the analysis period (dP 10% ). Therefore, we calculated the ratio between the mass-balance variations and increased precipitation (dB dP10% /dP 10% ) as an alternative sensitivity that indicates how much increased precipitation can enhance glacier mass gain. This ratio shows that a glacier receiving a moderate amount of precipitation (620-860 mm, Potanin and Turgen) will gain 1.6 to 2.0 times greater mass than that due to the increased precipitation. On the other hand, glaciers in arid environments (190-250 mm of annual precipitation, Sutai and Tsambagarav) have a much more sensitive response; 3.2 to 4.2 times greater mass gain than that due to the increased precipitation. The excess mass gain is caused by decreasing ablation with high albedo, which is kept by additional precipitation (Oerlemans and Fortuin, 1992;Fujita and Ageta, 2000). The contrasting responses suggest that the impact via albedo is highly affected by the wet or arid environments.  -295 ± 66 -248 ± 49 -292 ± 87 -496 ± 118 dB dP10% (mm w.e. dP 10% −1 ) 120 ± 17 61 ± 19 108 ± 36 138 ± 29 dP 10% (mm) 62 ± 6 1 9 ± 4 2 5 ± 4 8 6 ± 9 dB dP10% /dP 10% (-) 1.95 ± 0.32 3.19 ± 0.74 4.24 ± 1.55 1.64 ± 0.41 dB/dT σT (mm w.e.) -250 ± 56 -228 ± 45 -261 ± 77 -499 ± 119 dB/dP σP (mm w.e.) 112 ± 18 123 ± 28 166 ± 61 152 ± 38 a SMT is area-weighted summer mean temperature, P is annual precipitation, subscripts ERA5 and cal denote ERA5 and calibrated precipitation, B is glacier-wide mass balance, ELA is equilibrium line altitude, respectively. Significance level is p < 0.001 for all correlations (r) except for *p < 0.003, + p 0.2. dB/dT is mass balance sensitivity to temperature change. dB dP10% is mass balance sensitivity to 10% change in precipitation (dP 10% ). σT and σP are the interannual variability of annual mean temperature and annual precipitation, respectively.

Estimating the Mass Balance With Calibrated Precipitation
The precipitation calibration performed in this study is inevitably required owing to the limited observational coverage. However, it is likely that this calibration would still be necessary if there was a high degree of observational coverage, as the reanalysis product is unlikely to accurately represent accumulation-relevant processes, such as avalanching and wind scour/deposition. This is particularly evident across the Potanin and Turgen glaciers, where the stake locations were limited to the lower ablation zone, which led to the precipitation calibration largely affecting the glacierwide mass-balance estimation by altering the degree of accumulation at the higher elevations. Nakazawa et al. (2015) analyzed the pollen concentrations and oxygen stable isotopes in a snow pit on Potanin Glacier and determined the annual accumulation during the 2007-2011 period. The pollen-derived annual accumulations fit the linear regression lines of the stake-derived mass balance well, even though the timing of the mass balance measurements (June 2007-June 2008) was not consistent (Figure 3). However, a simple extrapolation of the linear regression would yield an unrealistically large amount of accumulation at the higher FIGURE 7 | Hypsometry and long-term mass balance profiles for Potanin (red), Sutai (yellow), Tsambagarav (blue) and Turgen (purple) glaciers. The stepwise lines denote the hypsometry at 50-m elevation bands. The thick solid lines denote the mean long-term profiles , and the associated shaded regions represent the ranges for each profile based on the tested precipitation parameters. The dashed lines denote the mass balance profiles resulted from precipitation gradients changing between 0 and 90% km −1 . Both amounts are larger than the calibrated annual accumulation values in this study (Table 3). We examined the water balances among the four glaciers (Supplementary Figure S9, annual mean glacier-wide values averaged for the 1980-2018 period). In contrast to the largely different water input (rain and snow, ∼600 mm w. e.), the difference of outgoing components (discharge and evaporation) falls within 300 mm w.e., which suggests that the contrasting glacier-wide mass balances are attributed to the amount of precipitation. On the other hand, the evaporations of Sutai and Tsambagarav Glaciers significantly contribute to the mass loss from the glaciers. This suggests that the sublimation estimation could be important for estimating mass balance and precipitation parameter of the glaciers under the fairly arid climate.
We validated the estimated mass balances with the remotely sensed surface elevation change rate (dh/dt, m a −1 ) provided by Hugonnet et al. (2021). The in-situ observed mass balances from Russian Altai and Tienshan were also compared (Supplementary Figure S10). We find that the decadal mean mass balances are not perfect but moderately consistent with the remotely sensed dh/dt if the density uncertainties (600-900 kg m −3 ) are taken into account. However, the data for Sutai are not consistent. We think that both data, our simulation and the dh/dt, could have large uncertainties. The short observation period (3 years) for estimating precipitation parameter of Sutai Glacier could cause an uncertainty. On the other hand, the mismatch for the west branch of Urumuqi No. 1 Glacier should be due to the uncertainty in the remotely sensed dh/dt, which would be mainly caused by weak surface contrast hampering generation FIGURE 9 | Climatic regime for the present-day geometry of (A) Potanin, (B) Sutai, (C) Tsambagarav and (D) Turgen glaciers for different summer mean air temperature (horizontal axes) and annual precipitation (vertical axes) combinations. Black open circles denote the long-term regime that was derived from the calibrated ERA5 data for the 1980-2018 period (39 years). Contour lines denote the glacier-wide mass balances in 0.5 m w. e. increments.
Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 785306 of digital elevation models. Different definitions and/or inconsistencies in the glacier outlines used in these studies could also affect the discrepancy. Supplementary Table S7 lists the water equivalent values of the estimated (B cal ) and observed (B obs ) mass balances against the remotely geodetic mass balance (B geo ). For converting the dh/dt value to water equivalent, we used the density of 850 kg m −3 which was assumed in Hugonnet et al. (2021). The biases (B geo -B cal or B obs ) suggest that the mass balance calculated in this study are estimated less negative than the geodetic mass balance (negative bias in B geo -B cal ) while the observed mass balance are more negative than the geodetic mass balance (positive bias in B geo -B obs ). Cross correlations support a plausibility of the reconstructed mass balance in this study (Supplementary We estimated the glacier-wide mass balances via linear regression (B reg ; blue lines in Figures 3-6), simulations using the uncalibrated ERA5 precipitation (B ERA5 ; red lines in Figures  3-6), and simulations using the calibrated precipitation (B cal ; black lines in Figures 3-6), and provided temporal change (Supplementary Figure  S11) and scatter plots (Supplementary Figure S12) for comparison. The estimated annual mass balance with a single precipitation parameter determined for each glacier is also depicted in Supplementary Figure S11 (B annual ). The annual mass balance is not always consistent with the calibrated mass balance (B cal ) because the precipitation parameter varies (Figures 3-6). The glacier-wide mass balances of the four glaciers with the calibrated precipitation exhibit negative mass balances throughout the in situ observation period (since 2004; thick lines in Supplementary Figure S11), whereas those based on a linear regression occasionally exhibit a positive mass balance (dashed lines). The scatter plots of glacierwide mass balance show that somewhat linear relationships approximately follow a 1:1 trend between the uncalibrated and calibrated precipitation results (orange circles in Supplementary Figure S12), whereas the linear-regression-based mass-balance results yielded rather scattered relationships (blue circles in Supplementary Figure S12). The comparison of these mass balances highlights inconsistencies in both the biases and trends of the calculated mass balances, which suggest that the glacier-wide mass balances of these four glaciers cannot be estimated without applying the calibrations performed in this study.
We assumed a precipitation gradient of 45% km −1 based on the elevation distribution of the precipitation parameter (Supplementary Figure S8), which is similar to a recently reported value for a Himalayan glacier (40 km −1 , Sunako et al., 2019) but much larger than those used for global simulations (Marzeion et al., 2012;Huss and Hock, 2015). On the other hand, precipitation could decrease with elevation where the vapor was not supplied (Putkonen, 2004). Therefore, the uncertainty in the precipitation gradient was tested by varying it over the 0-90% km −1 range (Figure 7). The reference elevation (z PR ) for the precipitation gradient was obtained from the mass-balance stakes, which were situated in the upper accumulation zone for the Sutai and Tsambagarav glaciers, and in the lower ablation zone for the Potanin and Turgen glaciers. Therefore, the change in the precipitation gradient affected the four glacier mass-balance profiles differently; a larger gradient yielded more negative mass balances for the Sutai and Tsambagarav glaciers, whereas it yielded more positive mass balances for the Potanin and Turgen glaciers (Figure 7). This variability in the precipitation gradient had a more significant impact on the Potanin and Turgen glacier mass-balance results, as more precipitation is expected on these glaciers than on Sutai and Tsambagarav (Table 3). However, the uncertainty due to the estimated precipitation parameter (r p ) is greater than that due to the gradient. In situ precipitation observations that cover a broad range of elevations are therefore required to determine the precipitation parameter and gradient. Furthermore, it is particularly important to measure the snowfall using weighingtype precipitation gauges.

Interannual Mass-Balance Variability
Significant correlations among the glacier-wide mass balance, SMT, and annual precipitation suggest that both parameters play key roles in the observed glacier mass-balance fluctuations. We multiplied the mass-balance sensitivity and interannual variabilities (standard deviation, σ) in SMT and annual precipitation to quantify the contributions of these two variables (Sunako et al., 2019). The variability in SMT exerts a greater effect on glacier mass balance than the variability in annual precipitation (dB/dT σT and dB/dP σP in Table 3, respectively), which range from 30 to 60% of the SMT influence. This difference could be attributed to the annual amount and/or interannual variability in precipitation. Furthermore, the temperature change could alter the melt amount in the ablation zone, while the precipitation change could alter the snow input in the accumulation zone. The bias induced by the changed variable could be further altered by the glacier hypsometry, which may vary considerably for top-heavy ice-cap and middle-heavy valley glaciers and would also alter the influence of these variables on the glacier-wide mass balance.

Equilibrium Conditions for the Glaciers
The zero isolines for the climate regimes and present-day hypsometries of the four glaciers (thick lines in Figure 9) indicate that the annual precipitation required to maintain equilibrium increases non-linearly with SMT, which is weighted by the hypsometry. The long-term means of the glacier-wide mass balances are negative for all of the glaciers ( Table 3 and Figure 10A).
The zero isolines appear along different sections of the four glaciers based on a range of SMT and annual precipitation conditions ( Figure 10A). We tested the zero-isoline sensitivity by exchanging the hypsometry and mass-balance profile to reveal the potential causes of these differences. The zero isolines of the valley-type glaciers occur close to each other when different glacier hypsometries are used with the Potanin Glacier massbalance profile ( Figure 10B). The hypsometry of the valley-type Turgen Glacier yields almost the same zero isoline, while the icecap-like Sutai and Tsambagarav glaciers yield slightly different curves. Significantly different zero isolines are obtained when different mass-balance profiles are applied to the Potanin Glacier hypsometry ( Figure 10C). These responses suggest that the climatic condition required for a given glacier to achieve steady state is not affected by the glacier geometry but rather significantly affected by the climatic condition itself.

Properties of the Mass-Balance Sensitivities
The gradient of the zero isoline is an indicator of glacier massbalance sensitivity (Ohmura et al., 1992;Braithwaite, 2008;Sakai and Fujita, 2017). Glaciers in cold, arid climates are generally less sensitive to temperature changes than those in warm, wet climates based on this gradient ( Figure 10A). The massbalance sensitivity to temperature changes follows this trend for the four studied glaciers (P cal and dB/dT in Table 3), whereby increased annual precipitation (P cal ) could result in a more negative mass balance for 1°C of warming (dB/dT). A comparison with similar numerical experiments (Oerlemans and Fortuin, 1992;Oerlemans, 2001) indicates that the mass balances of the four studied Mongolian glaciers are more sensitive to temperature changes than other glaciers worldwide. This can be attributed to seasonal precipitation variations, as summer accumulation is observed for the Mongolian glaciers while winter accumulation has been observed in previous studies (Fujita, 2008).
The mass-balance sensitivity to precipitation changes is more complicated than dB/dT ( Table 3). The abovementioned ratio (dB dP10% /dP 10% ) suggests that an increase in precipitation results in a more positive balance than the increased accumulation, especially in arid conditions (P cal and dB dP10% /dP 10% in Table 3). This increased precipitation could increase accumulation and cover the glacier surface with high-albedo snow that would suppress snowmelt. This process seems effective in cold, arid conditions. However, the increased precipitation in warm, wet conditions, such as those observed across Turgen Glacier, does not appear to induce an enhanced albedo effect; this is due to rainfall in a warmer environment. A comparison of these results with previous studies (Oerlemans and Fortuin, 1992;Oerlemans, 2001) suggests that the impact of precipitation changes across Mongolian glaciers is greater than those for glaciers in other regions; this could be also due to summer accumulation (Fujita, 2008). Significant precipitation changes occur during the summer, when both albedo and solid/ liquid phase changes alter the glacier mass-balance response.

CONCLUSION
This study reported the mass balances of four representative glaciers across the Altai Mountains in Mongolia, where in situ mass-balance observations have been obtained across the Potanin, Tsambagarav, Turgen, and Sutai glaciers since 2003, respectively. We estimated the precipitation parameters for the ERA5 precipitation values, and then evaluated the glacierwide mass balances using long-term mass-balance profiles that were calculated for the 1980-2018 period. Similar temporal fluctuations in the summer mean temperature anomaly and mass balance were observed, whereas the annual precipitation and ELA fluctuations were different among the glaciers. A crosscorrelation analysis suggested that these variables were strongly correlated with each other. We also tested the mass-balance sensitivity to temperature and precipitation changes. The sensitivities were correlated with the amount of annual precipitation, whereby wetter conditions resulted in a more negative mass balance for the same degree of warming, and an increase in precipitation had a larger impact on the mass balance in arid conditions. Continuous monitoring of both the glacier-wide FIGURE 10 | Balanced glacier-wide mass balances for Potanin (red), Sutai (yellow), Tsambagarav (blue) and Turgen (purple) glaciers. (A) Summary of Figure 9. Shaded regions denote the interannual variability of each zero-balance profile. Each point and its associated error bars represent the long-term mean and standard deviations of the calibrated ERA5 data shown in Figure 9, respectively. (B) Zero-balance isolines for the Potanin Glacier mass-balance profile with different hypsometries. (C) Zero-balance isolines for the individual mass balance profiles with the Potanin Glacier hypsometry. mass balance and local meteorological variables will contribute to improving our understanding of the state of mountain glaciers and their response to ongoing climate change across Mongolia. Furthermore, because the precipitation parameter and gradient could cause large uncertainty in the mass balance estimation, extending the stake network and precipitation observations via a weighing-type precipitation gauge would yield improved observations that possess lower uncertainties.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
PK designed and conducted the fieldworks. KF and AS designed the analysis of the study. PK and KF conducted mass balance simulations and analyzed the data. All authors contributed to discussion and writing the manuscript.

FUNDING
This study is supported by JSPS-KAKENHI (Grant No. 20H00196) and the Transnational Doctoral Program for Leading Professionals in Asian Countries in Nagoya University.