Adjustment of measurement errors to reconcile precipitation distribution in the high‐altitude Indus basin

Precipitation in the high‐altitude Indus basin governs its renewable water resources affecting water, energy and food securities. However, reliable estimates of precipitation climatology and associated hydrological implications are seriously constrained by the quality of observed data. As such, quantitative and spatio‐temporal distributions of precipitation estimated by previous studies in the study area are highly contrasting and uncertain. Generally, scarcity and biased distribution of observed data at the higher altitudes and measurement errors in precipitation observations are the primary causes of such uncertainties. In this study, we integrated precipitation data of 307 observatories with the net snow accumulations estimated through mass balance studies at 21 major glacier zones. Precipitation observations are adjusted for measurement errors using the guidelines and standard methods developed under the WMO's international precipitation measurement intercomparisons, while net snow accumulations are adjusted for ablation losses using standard ablation gradients. The results showed more significant increases in precipitation of individual stations located at higher altitudes during winter months, which are consistent with previous studies. Spatial interpolation of unadjusted precipitation observations and net snow accumulations at monthly scale indicated significant improvements in the quantitative and spatio‐temporal distribution of precipitation over the unadjusted case and previous studies. Adjustment of river flows revealed only a marginal contribution of net glacier mass balance to river flows. The adjusted precipitation estimates are more consistent with the corresponding adjusted river flows. The study recognized that the higher river flows than the corresponding precipitation estimates by the previous studies are mainly due to underestimated precipitation. The results can be useful for water balance studies and bias correction of gridded precipitation products for the study area.

Precipitation in the high-altitude Indus basin governs its renewable water resources affecting water, energy and food securities. However, reliable estimates of precipitation climatology and associated hydrological implications are seriously constrained by the quality of observed data. As such, quantitative and spatiotemporal distributions of precipitation estimated by previous studies in the study area are highly contrasting and uncertain. Generally, scarcity and biased distribution of observed data at the higher altitudes and measurement errors in precipitation observations are the primary causes of such uncertainties. In this study, we integrated precipitation data of 307 observatories with the net snow accumulations estimated through mass balance studies at 21 major glacier zones. Precipitation observations are adjusted for measurement errors using the guidelines and standard methods developed under the WMO's international precipitation measurement intercomparisons, while net snow accumulations are adjusted for ablation losses using standard ablation gradients. The results showed more significant increases in precipitation of individual stations located at higher altitudes during winter months, which are consistent with previous studies. Spatial interpolation of unadjusted precipitation observations and net snow accumulations at monthly scale indicated significant improvements in the quantitative and spatio-temporal distribution of precipitation over the unadjusted case and previous studies. Adjustment of river flows revealed only a marginal contribution of net glacier mass balance to river flows. The adjusted precipitation estimates are more consistent with the corresponding adjusted river flows. The study recognized that the higher river flows than the corresponding precipitation estimates by the previous studies are mainly due to underestimated precipitation. The results can be useful for water balance studies and bias correction of gridded precipitation products for the study area.

K E Y W O R D S
bias correction of precipitation, high-altitude Indus basin, net mass balance contribution to river run-off, net snow accumulation adjustments, precipitation distribution, precipitation measurement errors

| INTRODUCTION
High mountain ranges around the world are important sources of freshwater storage and subsequent supplies to downstream areas. Indus basin contains one of the most diversified and complex mountain terrains in the world. Precipitation in its high-altitude areas governs the renewable water resources determining water, energy and food securities in the region. Run-off regime of the basin is predominantly controlled by winter-and summer-monsoon precipitations and summer temperatures (Yu et al., 2013). Yet, there is limited understanding and reliable evidence of quantitative and spatio-temporal distribution of the key climatic variables, particularly the precipitation (Hewitt, 2005;Winiger et al., 2005;Ragettli and Pellicciotti, 2012;Immerzeel et al., 2015;Mishra, 2015) leading to a large uncertainty in the hydro-climatic predictability in the basin (Lutz et al., 2016). Overall scarcity and biased spatial and altitudinal distribution of the in situ observations are the primary reasons for this uncertainty and knowledge gap. Substantial increase in research on glacio-hydro-climatology of the Hindukush-Karakoram-Himalayan (HKH) region is observed since the International Panel on Climate Change (IPCC) released its fourth assessment report, which claimed that "glaciers in Himalayas are receding faster than in any other part of the world and, if the present rate continues, the likelihood of their disappearing by the year 2035 is very high" (Cruz et al., 2007). Later, IPCC withdrew this statement due to an inaccurate citation of the grey literature. Yet, most of the subsequent research is mainly focused on improved methods using more or less the same commonly available data sets that use low altitude and largely unrepresentative observations in the development or validation of these data sets.
Adequate monitoring of climatic variables to better represent the entire range of a diverse climate of this complex mountain terrain is essential for reducing uncertainties and inferring informed policy decisions. However, such an observational network in the study region is lacking mainly due to resource constraints and logistical limitations. To overcome the observational data gaps, the hydroclimatologists generally rely on numerous global/regionalscale gridded products derived through various means (e.g., climate models reanalysis, merged model and station observations, merged satellite estimates and station observations, and derived solely from station observations). However, the strong gradients and extreme heterogeneity of this complex mountain terrain are inadequately captured by the gridded products due to their coarse resolution and use of non-representative climate data in their development or validation (Immerzeel et al., 2015;Reggiani and Rientjes, 2015;Dahri et al., 2016). As such, the precipitation estimates by a number of earlier studies (e.g., Akhtar et al., 2008;Immerzeel et al., 2009;Bookhagen and Burbank, 2010;Bocchiola et al., 2011;Tahir et al., 2011;Immerzeel et al., 2012a;Mukhopadhyay, 2012;Central Water Commission and National Remote Sensing Centre, 2014;Lutz et al., 2014a;2014b;Reggiani and Rientjes, 2015) that used the gridded data sets show highly contrasting but consistently underestimated precipitation in most parts of the high-altitude Indus basin.
Numerous efforts to accurately estimate precipitation in this region only partially succeeded due to lack of observed data but significantly underlined the relevance and severity of the problem. In many hydrological modelling studies, the underestimated precipitation is often compensated for with other parameters like evapotranspiration and/or snow/glacier melt factors (Schaefli et al., 2005;Pellicciotti et al., 2012;Lutz et al., 2014a). This results in inaccurate and suboptimal inferences regarding precipitation distribution, snow/ glacier cover dynamics and associated melt water contributions. Adam et al. (2006) used a water balance approach to indirectly correct monthly precipitation in mountain regions from an existing global data set and provided reasonable approximations at basin level. However due to inaccuracies in water balance components and use of biased gridded data sets developed from limited observations, their results show large differences in precipitation amounts and distribution patterns at sub-basin scale in the study area. For example, precipitation in the high-mountain Karakorum region is largely underestimated due to lack of stations in this area, whereas higher precipitation amounts are shown for the southern parts of western Himalayan region that hosts many precipitation gauges. Lutz et al. (2014a) recognized underestimation of APHRODITE precipitation and multiplied it with an arbitrary constant factor of 1.17 to account for the inherent underestimations.
Recently, Immerzeel et al. (2015) and Dahri et al. (2016) used other sources of data/information to cover the observational gaps and provided relatively better estimates of precipitation amounts and distribution in the high-altitude Indus basin. The approach adopted by Immerzeel et al. (2015) used the glacier mass balance (GMB) estimates of Kääb et al. (2012) to inversely infer the high-altitude precipitation. Using APHRODITE as the basis, they computed vertical precipitation gradients until observed mass balance matched the simulated mass balance for the 550 major glacier systems in the Indus basin. However, precipitation in the basin does not have constant and linear gradients (Dahri et al., 2016), APH-RODITE precipitation distribution is highly biased (Palazzi et al., 2013;Dahri et al., 2016) and their mass balance computations are uncertain due to the use of extremely elusive direct evapotranspiration losses and negligence of percolation, interception and sublimation losses from the precipitation. Moreover, precipitation estimates of Immerzeel et al. (2015) might be affected by the overestimated basin boundaries of Shyok and Indus at Tarbela sub-basins. However, Dahri et al. (2016) integrated the available station observations with the indirect precipitation estimates at the accumulation zones of major glacier systems. They employed Kriging with external drift (KED) interpolation scheme with elevation as predictor to derive the spatio-temporal distribution of mean monthly and annual precipitation climatologies. They validated their precipitation estimates by the individual station observations and the observed specific run-off at sub-basin scale. However, if the net mass balance (i.e., slightly negative as estimated by Kääb et al., 2012) and precipitation losses (direct evapotranspiration, percolation, interception and sublimation) in the basin are taken into account, the Dahri et al. (2016) estimates still seem to be on lower side. The underestimated precipitation relative to the corresponding specific runoff in most sub-basins may be attributed to three possible reasons: (a) overestimated river flows, (b) significant contribution of snow/glacier melt without an adequate amount of precipitation to feed/sustain the glacier systems and (c) underestimated precipitation. Given the technological advancements and relative precision of discharge measurement techniques and quality control ensured by the data collecting agencies, river flows are generally considered to be adequately accurate. However, there is considerable speculation but little analysis and evidence regarding the contribution of net glacier mass imbalance to the river flows. Although Immerzeel et al. (2015) attributed the observed gap between precipitation and streamflow to the underestimated precipitation rather than the observed GMB, there is an emergent need to quantify the contribution of net glacier mass imbalance to the river flows. The underestimated precipitation by Dahri et al. (2016) is probably due to the use of net precipitation estimates from the glacier accumulation zones and the raw/uncorrected precipitation gauge observations which are subject to significant measurements errors (Sevruk and Hamon, 1984;Legates, 1987;Legates and Willmott, 1990;Goodison et al., 1998;Chen et al., 2015;Wolff et al., 2015).
The IPCC in its fifth assessment report stressed the need for adjustment of precipitation measurement errors and declared that observational uncertainties in precipitation may limit the confidence in the assessment of climatic change impacts (Bindoff et al., 2013). The measurement errors in precipitation observations, particularly the wind-induced undercatch of solid precipitation in windy conditions can be substantial (Adam and Lettenmaier, 2003;Wolff et al., 2015;Kochendorfer et al., 2017a;2017b). This is particularly important in the high-altitude Indus basin where moderately strong winds are a common phenomenon; temperature mostly remains below the freezing point and the majority of precipitation falls in the form of snow. Legates (1987), Legates and Willmott (1990) and Adam and Lettenmaier (2003) adjusted the systematic biases of global precipitation products including the Indus basin but these data sets included only a few stations located in relatively dry valleys in the study area. The uncertainties in precipitation estimates may significantly affect the outcomes of hydrological/land surface models and mass balance studies. A systematic error of over 3% in rainfall measurement could lead to substantial underestimation of water in the hydrologic system (e.g., Sevruk, 1982;Biemans et al., 2009). Therefore, the systematic errors in precipitation observations must be corrected if the measurements are to be used for climate change, hydrological modelling and water balance studies (Legates and Willmott, 1990;Voisin et al., 2008;Wolff et al., 2015). This study attempts to address the above concerns by adjustment of the systematic measurement errors in precipitation observations, adjustment of net snow accumulation for the ablation losses and adjustment of river flows for the net mass balance contributions. The ultimate goal of this research is to facilitate creation of an accurate and consistent gridded precipitation product for the highly under-explored region of Indus basin. The results will have considerable implications for water resources planning and management in both upstream (high altitude) and downstream (low altitude) areas of the Indus basin.

| STUDY AREA
The study area covers the high-altitude catchments of the Indus river, which originates from the Tibetan Plateau (TP) and the HKH mountain regions ( Figure 1). The total area of the study region is about 4.03 × 10 5 km 2 of which 50% is above 4,000 m a.s.l. and another 24% between 2,500 and 4,000 m a.s.l. Precipitation in the study area is influenced by multiple weather systems. The Indian summer monsoon brings moisture from the Indian Ocean and Bay of Bengal and is the dominant system in the southeastern areas. The western disturbances originating from the Mediterranean and Caspian Sea dominate the southwestern and northwestern areas bringing winter monsoon during December-April months. During spring and early summer, irregular collapses of the Tibetan anticyclone sometimes allow monsoonal air masses to penetrate into the Karakoram Range (Wake, 1989). Direct transport of moisture from the Arabian Sea and local evapotranspiration also have considerable influence as about 5-40% of the precipitation falling in the Himalayas originates from the irrigated areas in northern India and Pakistan (Tuinenburg et al., 2012;Harding et al., 2013;Wei et al., 2013). However, the hydrological cycle in the study region is usually intensified when all or some of these systems interact with each other.
Information regarding the gauge type, use of wind shield if any, orifice area and height of the gauge orifice were taken from Sevruk and Klemm (1989) and Bureau of Indian Standards (1992a, 1992b and from PMD and WAPDA through personal communications. Until 1969, the most extensively used rain gauge in India was nonrecording (Symon's gauge or MK2 model) with orifice area of 127 cm 2 and instrument height of 0.3 m (Sevruk and Klemm, 1989). Thereafter, Indian standards adopted by the Bureau of Indian Standards (BIS) for design and manufacturing of meteorological instruments are strictly followed and Indian rain gauge (20-22-P) reinforced with fibreglass polyester is predominantly used (Bureau of Indian Standards, 1992a;1992b). Similarly, the most widely used rain gauge type by PMD has been nonrecording MK2 (13-15-C) model with orifice area of 127 cm 2 and instrument height of 0.3 m. In 2010, PMD started using its own model, which is tipping bucket rain gauge (TBRG) type equipped with logger and standalone method of monitoring rainfall, with 0.2 mm (moderate rain) tipping bucket, orifice area of 400 cm 2 and gauge height of 0.6 m. WAPDA uses both automatic weighing and standard meteorological service manual rain gauges. The automatic gauges have an orifice area of 127 cm 2 , tipping capacity of 0.254 mm and gauge height of 0.3 m (Water and Power Development Authority, 2008). A manual gauge is read in conjunction with each automatic gauge as a check on the total rainfall. In 1994-95, WAPDA installed 20 automatic data collection platforms (DCPs) in the high-altitude areas that use snow pillows to measure both solid and liquid precipitation as water equivalent (SIHP, 1997). The observatories installed and maintained by the University of Bonn under the CAK program used the automatic weather stations including data logger, tipping bucket and snow depth gauge to measure precipitation (Miehe et al., 1996). Afghanistan mainly uses the Tretyakov (20-24-G) type of rain gauge without windshield having orifice area of 200 cm 2 and 0.4 m height (Sevruk and Klemm, 1989). The metadata of 305 precipitation observatories and 21 glacier observation points used in this study are outlined and described in Table S1, Supporting information).

| Temperature and wind speed observations
The adjustments for wind-induced under-catch of precipitation observations require corresponding data of temperature and wind speed. However, out of 324 stations, temperature data were available for only 114 stations (Table S1). We therefore derived monthly lapse rates based on elevation and latitude and estimated the maximum and minimum temperatures for the remaining stations. The observed data of wind speed was available for only 25 stations. Wind speed for the remaining stations is taken from the Japanese 55year Reanalysis (JRA55) data set (Kobayashi et al., 2015). JRA55 provides wind speed estimates at the standard anemometer height of 10 m, whereas the station-based observed wind speed is measured at 2 m height. In order to get an idea of the accuracy of the JRA55 wind speed data, we compared it with the observed wind speed for the 25 stations. For this purpose, we computed wind speed from the U-and V-components at 10 m height and downscaled it to match the 2 m height of stations using the Monin Obukhov theory (Businger and Yaglom, 1971;Obukhov, 1971). Although we could not detect large differences and/or any definite and strong trends, a tendency of slightly underestimated wind speed in low-altitude areas and vice versa in high-altitude areas is noticed. We also observed marginally increased wind speeds during November-February months and slightly decreased wind speeds during March-October months for the JRA55 data. Due to insufficient observed data of wind speed, we have neglected these minor differences and used wind speed data of JRA55 as such. Nevertheless, such minor differences of wind speeds in JRA55 data might result in slight overestimation of precipitation adjustments in the higher-altitude areas during four (November-February) winter months and slight underestimation of precipitation adjustments in the lower-altitude areas during the remaining months.

| River flows
Daily data of the observed river flows at sub-basin level for 14 hydrological stations ( Figure 1) in the study area were collected from WAPDA. We used flow data of Jhelum and Chenab rivers for 1961-1970 period and all the rivers in the western part sub-basins for 1999-2011 period to coincide with the precipitation data periods. Ravi, Beas and Sutlej basins are located in India and their inflow data are not publicly available. Therefore, we extracted mean monthly river flows from Adeloye et al. (2016)  It is worth to note that there are considerable diversions in some sub-basins on the upstream side of their rim stations (e.g., at Warsak, Nowshera and Tarbela), which are often overlooked by previous studies. We also collected the data of these upstream diversions and added them to the flows of the respective sub-basins. River flow data of coinciding time periods are used to validate the adjusted precipitation at sub-basin scale.

| Precipitation measurement error adjustment methods
The amount of actual precipitation reaching the ground is generally higher than what is measured in precipitation gauges due to measurement errors, which usually depend on the form of precipitation, gauge type, topography, vegetation around the gauge site and the exposure of the gauges to prevailing temperatures and winds. Wind-induced undercatch is by far the most dominant source of errors in gaugemeasured precipitation observations (Goodison et al., 1998;Adam and Lettenmaier, 2003;Michelson, 2004;Wolff et al., 2015), yet most of the widely used global precipitation data sets are not adjusted for such errors (Adam and Lettenmaier, 2003). While recognizing the significance of measurement errors in precipitation observations, the World Meteorological Organization (WMO) initiated a comprehensive program of international precipitation measurement intercomparisons during 1960-1993 and established the pit gauge (Sevruk and Hamon, 1984) and the double-fence international reference (DFIR) (Goodison et al., 1998) as the standard reference gauges for liquid (rain) and solid (snow) precipitation, respectively. Sevruk and Hamon (1984) and Goodison et al. (1998) also underlined the need for gauge calibration and adjustment of errors to increase reliability of the precipitation data. However, the agencies involved in measurement of precipitation in the Indus basin generally indicate to follow the WMO standards for design, construction, installation and operation of precipitation gauges but hardly or inadequately adjust the systematic measurement errors at the source, which signifies the need for correction of measurement errors. Sevruk (1982) related and statistically analysed various components of the systematic measurement errors to the meteorological and instrumental factors and proposed a general equation for adjustment of gauge-measured precipitation errors. Legates (1987) later modified it to account for both liquid and solid precipitation components separately. The modified equation is expressed as P a = 1−R ð ÞK r P m +ΔP wr +ΔP tr +ΔP er ð Þ +RK s P m +ΔP ws ð where P a is adjusted precipitation (mm), R is proportion of solid precipitation, K is correction coefficient that accounts for wind-induced losses, P m is measured precipitation (mm), ΔP w is wetting losses (mm), ΔP e is evaporation losses (mm), ΔP t is trace precipitation (mm) and subscripts r and s denote rain and snow components, respectively. Legates (1987) model was developed for a variety of manual rain gauges including Nipher, Tretyakov and MK1/MK2 models with and without windshields. However, significant uncertainties remained for wind-induced under-catch of solid precipitation particularly by automatic precipitation gauges. Nitu and Wong (2010) observed much larger variation between gauges and windshield configurations for automatic stations than for manual stations. Wolff et al. (2015) compared precipitation data from the standard automatic Geonor precipitation gauge with data from a reference configuration consisting of an automatic precipitation gauge (Geonor T200-BM) and an Alter wind shield with double-fence construction. They derived an adjustment model to determine catch efficiency as a continuous function of both wind speed and air temperature using Bayesian statistics to more objectively choose the model that best describes the data. Wolff's model allows solid precipitation adjustments at wind speeds greater than 7.0 m/s. However, it is also gauge/shield-specific and different site specificities and gauge/shield configurations might result in different adjustment functions. Kochendorfer et al. (2017a) analysed precipitation measurements from eight different WMO-SPICE sites for both unshielded and single-Alter-shielded OTT Pluvio 2 and Geonor T-200B3 types of weighing gauges. They grouped unshielded and single-Alter-shielded precipitation gauge configurations separately irrespective of gauge types and created a single transfer function of air temperature and wind speed using the corresponding measurements from the reference gauge. They also derived the coefficient fits for both unshielded and single-Alter-shielded precipitation gauges at gauge height as well as 10 m height. The derived transfer function is expressed as where T air is mean air temperature ( C), U is wind speed (m/s), a, b and c are the coefficients fit to the data and TAN −1 is the inverse of tangent function. Our method of adjusting systematic errors in precipitation measurements largely follows the approach by Adam and Lettenmaier (2003) using the "liquid" part of the model by Legates (1987) and uses the model by Kochendorfer et al. (2017a) for adjustment of the solid precipitation component. The detailed methods for computation of the required variables in Equation (1) are described in the supplement available on-line. The coefficient values in Equation (2) (a = 0.0623, b = 0.776, c = 0.431) are taken as determined at 10 m height by Kochendorfer et al. (2017a). We used the coefficient values of 10 m height because most of our wind speed data belonged to the JRA55 data set, which provides wind speed data at 10 m height. The observed wind speed at 25 stations is converted from observation height to 10 m height using the Monin Obukhov theory (Businger and Yaglom, 1971;Obukhov, 1971).

| Adjustment of net snow accumulations
The meteorological stations in the study area are unevenly distributed in both horizontal and vertical direction. Scarcity of precipitation measurements at higher-altitude areas, where the bulk of precipitation falls, seriously limits an accurate assessment of precipitation climatology and its hydrological implications. In order to overcome this observational data gap, we assumed 21 virtual stations at the major glaciers where the net snow accumulations were estimated through mass balance studies using snow pillows, snow pits and ice cores (e.g., Qazi, 1973;Decheng, 1978;Batura Investigations Group, 1979;Kick, 1980;Mayewski et al., 1983;Mayewski et al., 1984;Wake, 1989;Bhutiyani, 1999;Shroder et al., 2000;Mayer et al., 2006;Hewitt, 2011;Mayer et al., 2014). However, most of these mass balance studies were undertaken in the active ablation zones of the glaciers, where ablation and accumulation processes are happening simultaneously. Generally, glacier ablation is the function of ablation rate, altitude of the equilibrium line altitude (ELA) and the elevation difference between mean ELA and the glacier observation point. Ablation zones are the areas below the ELA, which is the elevation at which the annual net mass of the glacier remains zero and the area above this elevation is known as the accumulation zone (Cuffey and Paterson, 2010). Hence, the estimated net glacier mass accumulations are subject to ablation losses until the next accumulation period. The ablation gradients can be variable depending on debris cover and surface albedo or energy availability to melt the exposed glaciers. Wagnon et al. (2007) observed ablation gradients of 0.60-0.81 m w.e. (water equivalent) for each 100 m with a mean value of 0.69 m w.e. over a period of 4 years of mass balance studies at the Chhota Shigri glacier, western Himalaya. Yu et al. (2013), based on glacier studies by Mayer et al. (2006) and Wagnon et al. (2007) in the Karakoram and western Himalaya, assumed an ablation gradient of 1 m w.e. per 100 m for the upper Indus basin. Hewitt et al. (1989) however, estimated an ablation gradient of 0.5 m per 100 m for the middle portion of the ablation zone on the Biafo glacier in the central part of the Karakoram. No ablation above ELA is assumed. We selected the rather conservative estimates of ablation gradient by Hewitt et al. (1989) and adjusted the net accumulations by taking the ELA as the boundary for the ablation process. However, the location of ELA can vary from location to location. In temperate glaciers, usually the snow line elevation (SLE) and ELA are often assumed to be the same. The estimates for mean ELA at sub-basin scale are taken from Khan et al. (2015), who estimated ELA values based on SLE.
3.6 | River flow adjustments WAPDA uses standard flow measuring devices to ensure high quality river flow data. The primary river flow measuring technique uses area velocity measurements to determine the stage-discharge relationships and associated rating tables. The results are verified by area-velocity method, area-slope method, contracted opening measurements, or computation of flow over dams or weirs (Water and Power Development Authority, 2012). The daily mean discharge values are computed from the mean gauge heights and corresponding calibrated rating tables. In case of extremely high discharges, the rating curves are extrapolated by applying simple linear regression between the gauge height and discharge measurements. The actual measurements are however taken 4-8 times per month. The intermediate daily values are estimated from the rating tables. The accuracy of stream flow measurements depends primarily on stability of the stage-discharge relationship, frequency of discharge measurements if the relationship is unstable, and accuracy in the observation of the stage and measurement of discharges. In general, monthly and annual mean values are more accurate than daily values because of compensation of random errors. WAPDA evaluates the probable accuracy of discharge measurements as excellent (error < 5%), good (error < 10%), fair (error < 15%) and poor (error > 15%). In general, a probable accuracy of 0-5% is aimed for.
Although river flow data may still be subject to some degree of uncertainty due to measurement errors, we assumed river flows as adequately accurate considering the relative precision of discharge measurement techniques and quality control ensured by the data collection agencies.
To account for the contribution of net glacier mass imbalance in each sub-hydrological basin, we adjusted the measured river flows. Kääb et al. (2012) used satellite laser altimetry and a global elevation model and observed a slightly negative mass balance of −0.21 ± 0.05 m/year w.e. for HKH region during 2003-2008 with maximum rates of −0.66 ± 0.09 m/year w.e. in the western Himalayan (Jammu-Kashmir) areas. We derived the specific net mass balance rates at sub-basin scale from the mass balance estimates of Kääb et al. (2012) and took glacier areas from the Randolf Glacier Inventory (RGI) version 5.0 (Arendt et al., 2015) to compute the contribution of the changes in the net glacial mass imbalance to the observed river flows. The adjusted river flows are used for validation of the adjusted precipitation estimates at sub-basin scale.

| Spatial interpolation
The actual and error-adjusted point measurements of mean monthly precipitation are spatially interpolated following Dahri et al. (2016), who used the KED interpolation scheme (Schabenberger and Gotway, 2005) with elevation as a predictor to derive spatio-temporal distribution of precipitation in the high-altitude Indus basin. The KED model includes a component of spatial autocorrelation and a component for multi-linear dependence on pre-defined variables (predictors). It considers the observations (Y) at sample locations (s) as a random variable of the form (e.g., Diggle and Ribeiro, 2007), where μ(s) describes the deterministic component of the model (external drift or trend) and is given as a linear combination of K predictor fields x k (s) (trend variables) plus an intercept (β 0 ). The β k are denoted as trend coefficients, while Z(s) describes the stochastic part of the KED model and represents a random Gaussian field with a zero mean and a second-order stationary covariance structure. The latter is conveniently modelled by an eligible parametric semivariogram function describing the dependence of semivariance as a function of lag (possibly with a directional dependence). Dahri et al. (2016) provided a detailed account of the KED interpolation method including model description and functionalities, reasons for its selection and comparative advantages of its use in the high-altitude Indus basin.

| Cross validation of the adjusted precipitation
We used exactly the same approach of interpolation and cross validation as adopted by Dahri et al. (2016), where the cross validation applied on the observed and predicted values from all the stations is used to assess the errors/ uncertainty associated with the interpolation scheme by using error scores of the relative bias (B) and the relative mean root-transformed error (E), which are defined as where P i and O i are the predicted and observed precipitation values, respectively, while O is the average of all (or a subset of ) the station observations and n refers to the number of precipitation values. Under ideal conditions, the overall performance of the employed regression models and interpolation estimates at basin/sub-basin scale can also be cross validated by applying the continuity equation suggested by Budyko (1974), which is given by where P, Q, ET and G are the basin-average precipitation, run-off, evapotranspiration and net groundwater discharge, respectively, while ΔS is the net change in storage for a given time increment (Δt). Equation (7) can be modified by adding interception (I), sublimation (S) and net mass balance (ΔMB) contributions for the highly glacierized and snowpack-dependent river basins as follows: Unfortunately, there are no independent data sets of actual evapotranspiration, sublimation, interception and the net groundwater discharge for the study area. The globalscale data sets of these variables are generally more uncertain than precipitation itself; therefore, it would be unwise to validate the estimated precipitation with these extremely uncertain data sets. Nevertheless, surface storage and groundwater recharge are mostly very low in high-altitude areas, which are mostly rocky bare mountains with steep slopes and no groundwater. Precipitation may travel long distances through breaches but ultimately joins the river streams as base flow. Although there might be considerable delay effects, these may be considered negligible for longterm average conditions. Similarly, the surface storage due to topographical undulations may also have a delaying effect. Interception by the vegetation cover and sublimation (direct evaporation from the snow glacier fields) are included in the total direct evapotranspiration. Direct evapotranspiration is notoriously complex to measure as it is among others a function of water availability as well as water demand. The available global-scale gridded data sets of actual evapotranspiration are highly inconsistent in quantitative as well as spatial distribution terms and generally reflect overestimated values. We therefore rely mainly on the specific run-off and net mass balance data to validate our adjusted precipitation estimates.

| Precipitation adjustments
To facilitate adjustment of measurement errors in precipitation observations, the corresponding air temperature is determined from elevation and latitude based lapse rates. The results revealed a strong correlation of temperature with elevation and considerable correlation with latitude ( Figures S2-S5). Significantly different gradients for each month and substantial difference among the gradients for maximum and minimum temperatures were observed (Table 1). Hence, use of a universally assumed or timeindependent site-specific observed gradient of mean annual temperature to estimate maximum and minimum temperatures (e.g., Immerzeel et al., 2012a; 2012b; Lutz et al., 2013) is probably not correct in the high-altitude Indus basin. Comparison of Table 1 and Figures S2 and S3 indicate that incorporation of latitude as an additional predictor improves the correlation of the regression models by up to 6.0% for maximum temperature and up to 1.5% for minimum temperature during 1999-2011. Almost similar trends are observed for 1961-1970 period. The contribution of elevation to the correction is positive in the summer months and negative in the winter months, while the contribution of latitude is positive throughout the year. The highest improvement is achieved during the monsoon season (July-September).
To illustrate the precipitation biases over the highaltitude Indus basin, the results for each individual station are presented. The applied bias adjustments significantly increased the gauge-measured precipitation. The highest increments are computed for wind-induced under-catch of solid precipitation followed by liquid precipitation undercatch, wetting losses and precipitation losses during trace events (Figure 2a-d). The solid precipitation under-catch generally dominates the higher-altitude stations, that is, elevations greater than 2000 m and during the December-April months. The range of liquid precipitation undercatch is much lower and mainly concentrates in the summer monsoon dominated low-altitude areas, that is, elevation less than 3,500 m. The wetting losses and unmeasured trace precipitation depend on the number of precipitation events. In many cases, particularly for the low-altitude stations experiencing lower wind speeds, the wetting losses exceeded the wind-induced under-catch of liquid precipitation due to the fact that it covers all the stations and both forms of precipitation (liquid and solid). The total bias between the gauge-measured and erroradjusted precipitation ranged from 12 to 773 mm/year for various individual stations and up to 1,000 mm/year for the glacier points ( Figure 2e). The total absolute biases (corrections) for all the stations at monthly and annual scale are given at Table S2. The largest increases are found for the stations receiving greater precipitation amounts, located at higher-altitudes and encountering higher wind speeds. Based on the above mentioned corrections, we introduced monthly-scale correction factors (CFs) for each station (Table S3). These station based CFs vary over space and time, with stronger magnitude in higher-altitude areas ( Figure 2f ) and during winter months (Table S3).

| Snow accumulation adjustments
The total ablation losses at a given ablation rate from a glacier zone depend on the ablation gradient and ΔELA (the difference between the mean elevation of a glacier zone and ELA). Assuming that the practical ablation above ELA is insignificant, the potential ablation losses from the selected glacier zones vary from 0 to 1,000 mm/year (Table 2). These ablation losses are added to the original estimates of the net accumulations to account for the ablation losses from the actual precipitation.

| Spatial distribution of unadjusted and adjusted precipitation
Continuous fields of precipitation generated through KEDbased interpolation of the adjusted station observations and adjusted snow accumulations at monthly scale show how precipitation patterns and amounts are spatially distributed in the study area (Figure 3a-l). Monthly precipitation distributions largely confirm the bimodal weather system reflecting the wintertime precipitation associated with the westerlies and the impact of Indian summer monsoon in the study area. Overall climatology and distribution patterns of the adjusted precipitation (Figure 3m) match very well to the unadjusted case (Figure 3n) or estimates of Dahri et al. (2016). However, the adjustments revealed significant improvement in terms of quantitative and spatio-temporal distribution of precipitation in the study area (Figure 3o). An overall increase of 21.3% in average annual precipitation is realized at basin (study area) level, while at sub-basin scale it ranged from 6 to 77% (Table 3). Greatest improvements are achieved in the high-altitude areas of Astore, Shyok, Shigar, Hunza, Gilgit and Chitral sub-basin and during the winter months.

| River run-off adjustments
The net mass balance estimates of Kääb et al. (2012) for the study area are translated into the amount of run-off generated at sub-basin scale. As a result of slightly negative mass balance estimates for all sub-basins, their contributions to river run-off are also negative and relatively small ranging from 0.4 to 6.1%. The adjustments in river-specific run-off depend on the net mass balance as well as glacier area and varied from −51.5 mm in the Chenab sub-basin to −2.5 mm in the Panjkora sub-basin (Table 4).

| Validation of precipitation estimates
The estimated precipitation distributions can be validated by evaluating the accuracy of the employed interpolation scheme and the output interpolated fields. For accuracy assessment of the interpolation scheme, the KED interpolation model produces both prediction as well as error/uncertainty surfaces, giving an indication or measure of how good the predictions are. The cross validation applied on the observed and predicted values from all the stations resulted in relative bias (B) error scores of less than 1, suggesting a negligible underestimation of the predicted values for all months except August, which shows a slight overestimation (Table 5). Similarly, the relative mean roottransformed error (E) scores of less than 1 for the months January-May suggest excellent results. While the remaining months of June-December experience E values of greater than 1, which depict typical errors slightly greater than the spatial variations. Almost similar trends are observed for the unadjusted case. In general, the cross-validation results  Figure S1) covering the two time periods of 1999-2011 and 1961-1970, respectively. T x 1-12 and T n 1-12 refer to the calendar months for maximum and minimum temperatures, respectively. E denotes elevation (m) and L represents latitude (decimal degrees) of the meteorological stations. R 2 is the combined correlation of temperature with E and L Another means of validation is the comparison of the estimated precipitation with the corresponding observed river flows (specific run-offs). Dahri et al. (2016) demonstrated that the previous estimates of precipitation distribution in the study area are not only highly contrasting but largely underestimating the actual precipitation. Likewise in the Dahri et al. (2016) study, precipitation estimates derived from the unadjusted precipitation observations provided relatively better estimates than the previous studies. Yet, slightly lower precipitation than the measured specific runoff in 9 out of 17 sub-basins (Figure 4) is absolutely counterintuitive implying underestimated precipitation or an unaccounted source of water (e.g., glacier melt contribution). Long-term annual mean precipitation must always be greater than the corresponding specific run-off if a positive or neutral mass balance is prevalent in any river basin. In case of a negative mass balance, its contribution to river flows has to be subtracted from the actually observed river flows and the adjusted flows must be lower than the corresponding mean annual precipitation. Cross validation of adjusted precipitation estimates with the corresponding adjusted specific run-offs (Figure 4) revealed adjusted specific run-off well below the adjusted precipitation estimates for all the sub-basins except Swat, which reflects underestimated precipitation or a bigger contribution of a negative mass balance to river flows.

| DISCUSSION
Precipitation is an integral component of the hydrological cycle and usually the most important input to water balance assessments and climate change studies. Hence, its accuracy is essential as errors in precipitation estimates may translate into major changes in the water budget of a particular region. However in many areas, precipitation measurements are still subject to significant errors and a large uncertainty  05 1.1 1.15 1.2 1.3 1.4 1.5 1.8 2.0 2.3 2.6 3.0 precipitation. The situation is particularly serious in the high-altitude Indus basin where biased distribution and lack of the observed data further worsen the problem. As such the precipitation products derived from or validated by the observed data covering this region are prone to significant errors (Reggiani and Rientjes, 2015;Dahri et al., 2016).
Scientists have used different approaches to overcome the observational data gaps. For example Adam et al. (2006) used a water balance approach to indirectly estimate precipitation. However, large uncertainties in the different water balance components limit wider application of this approach. Immerzeel et al. (2015) used mass balance estimates to inversely compute precipitation in the major snow/ glacier zones and applied a linear lapse rate of precipitation increase with elevation up to 5,000 m using APHRODITE as the reference data set. Uncertainties in mass balance and water balance components and assumption of linear precipitation increase with altitude are the major drawbacks of this method. Dahri et al. (2016) integrated station observations with the net snow accumulations estimated through mass balance studies and applied KED interpolation scheme to derive precipitation in ungauged areas. Measurement errors in station observations and negligence of snow/glacier ablations in the net snow accumulations are the key shortcomings of this approach.
The approach adopted in this study is based on catch adjustments of precipitation observations for systematic measurement errors, adjustment of net snow accumulations for the ablation losses and adjustment of river flows for the contribution of net GMB. Mean monthly precipitation climatologies are derived from the actual precipitation observations and actual net snow accumulations as well as from the adjusted precipitation observations and the adjusted net snow accumulations following Dahri et al. (2016).
The results presented in this study further support the wind-induced under-catch as the largest source of errors in gauge-measured precipitation observations. The catch corrections have increased the gauge-measured precipitation values ranging from 12 to 773 mm/year for various stations, while net snow accumulations at the glacier points increased up to 1,000 mm/year. A large part of precipitation in the high-altitude Indus basin falls as snow, which is more susceptible to under-catch even at moderate wind speeds. The largest corrections were found for wind-induced under-catch of solid precipitation, which is in line with the results of previous studies (e.g., Legates and Willmott, 1990;Goodison et al., 1998;Adam and Lettenmaier, 2003;Michelson, 2004  2 Adjusted net snow water equivalent at the major glacier accumulation zones. Lon. is longitude, Lat. is latitude, Ele. is elevation, ELA is equilibrium line altitude, ΔELA is the net elevation contributing to ablation and ΔA is adjustment in the net accumulation loss and trace precipitation loss are also important, particularly in low-altitude and relatively dry areas. The large differences between the observed precipitation and the corresponding specific run-off observations (usually greater specific run-off than precipitation) in previous estimates are often attributed to the contribution of snow/glacier melt. Indeed the high-altitude Indus basin receives considerable snow/glacier melt contributions, which largely come from the melting of temporary/seasonal snow cover and may vary from year to year depending on the quantity and timing of winter snowfall and snowmelt during the subsequent summer. However, quantitative estimates of net GMB contributions to river flows are largely lacking. Therefore, the accuracy of the estimated net GMB contributions to the river flows is mainly depending on the uncertainties in glacier area and the ablation rates of mass balance. Our methodology of adjusting river flows for the net mass balance contributions is straight forwards and the adjustments are slightly less than what is modelled by Lutz et al. (2016). For example, we estimated net GMB contribution of  Lutz et al. (2016). The difference might be due to the use of different approaches and different glacier inventories having different glacier areas. Lutz et al. (2016) pointed out a 23% difference in the glacier areas from three different inventories implying considerable differences in the water balance components.
The precipitation distribution derived through actual station observations combined with the actual net glacier accumulations is almost similar to that derived by Dahri et al. (2016) except for the addition of a few sub-basins and the use of additional and updated observed data. The catch corrections and snow accumulation adjustments significantly increased the total gauge-measured as well as basin-scale precipitation (Figures 2-3o and 4 and Table 3). The overall distribution patterns of precipitation remained largely the same as identified by Dahri et al. (2016), but substantial increases in the magnitude of precipitation amounts are realized. One of the advantages of the KED interpolation method is that it estimates an interpolated surface from a randomly varied small set of measured points and recalculates estimated values for these measured points to validate the estimates and determine the extent of errors. When compared with the corrected precipitation derived by Immerzeel et al. (2015), our estimates show significantly smaller root-mean-square error and a stronger correlation with the error-adjusted station observations ( Figure 5). The corrected precipitation estimates by Immerzeel et al. (2015) show considerable differences with significantly lower values at the majority of station locations including the points at the major glaciers, where actual measurements of net snow accumulations were taken. At the basin scale their estimates are relatively better but seem to be on the higher side in about half of the sub-basins. This discrepancy between station-based point observations and basin-scale precipitation estimates by Immerzeel et al. (2015) may be attributed to the higher and linear lapse rates of precipitation increase applied to compute the precipitation fields. Also, they did not validate their estimates with the observed precipitation of the individual stations. Instead, they used the Turc-Budyko representation to show the physical realism of their estimates and attributed some of the estimates that fall on the right side (inside) of the theoretical Budyko curve to the possible contribution of the negative mass balance to river flows and uncertainties in the potential evapotranspiration (ET p ) data set.
In this study, we used accurate run-off observations (specific run-offs), which are further improved by adjusting for the net GMB contributions, and improved ET P estimates from JRA55 reanalysis data set ( Figure 6) to evaluate the physical realism of our estimated precipitation compared to the precipitation estimates from Immerzeel et al. (2015). Over one third of the points representing estimated precipitation by Immerzeel et al. (2015) in various sub-basins (e.g., Gilgit, Chitral, Panjkora, Kabul at Warsak and Nowshera, and Sutlej) lay inside the theoretical Budyko curve indicating higher values than the theoretically expected. However, the estimates of unadjusted precipitation in our study, which are almost similar to the estimates of Dahri et al. (2016), show 10 out of 17 sub-basins above the line of moisture limit indicating underestimated precipitation in these sub-basins. The adjusted precipitation derived in our study shows relatively better fits in the Turc-Budyko representation except for the Swat sub-basin. The greater specific run-off than precipitation in the Swat basin may be attributed to yet an underestimated precipitation and/or greater negative mass balance than what is presently assumed.
The run-off ratio (Q/P) determines the amount of precipitation converted into overland flow or surface run-off. It is mainly controlled by largely stable natural factors including climate, soil and topography and to some extent by the human alterations to landscapes. Relatively higher run-off ratios are produced for areas with shallow or clay soils, steeper slopes and devoid of vegetation cover. Snowcovered areas hold winter precipitation as snow/ice and  Figure S6). All these topographical properties infer the high-altitude Indus basin as a typical case of an area that accelerates rapid run-off generation. Therefore, relatively high rates of run-off ratios are to be expected. Table S4 and Figure 6 show the improved run-off ratios (Q/P) and aridity indices (P/ET p ) if compared to the data sets of Immerzeel et al. (2015) and Dahri et al. (2016). Although the error-adjusted precipitation derived in this study seems to be more consistent, yet there are a few uncertainties that need to be understood and taken care of in future investigations. The major uncertainties associated with the results of our study may arise from four possible sources: (a) uncertainties in regression models due to their imprecision and uncertainties in the input data, (b) uncertainties arising from the estimated temperature and wind speed for many observatories, (c) uncertainty in the gauge type of the basin's gauge network and (d) uncertainties in spatial interpolation of the point observations to derive gridded fields of precipitation. The error estimation of the regression models employed in this study are tested at different locations and the relationships with the best fit are also applicable for similar situations in other areas. Nevertheless, regression models are in essence approximations of reality and some degree of uncertainty will always remain in the results. Relatively more accurate adjustments of precipitation under-catch for any precipitation event can be made by using the corresponding data of temperature and wind speed. However, hourly or daily data of these parameters are not available for many observatories in the study area. Also, there are many stations for which such data are not available at all. For locations without these data, temperature may be derived from the lapse rates of the available observations and wind speed from JRA55 data set. However as shown, the use of these data may add to the uncertainties in the catch corrections. The meteorological data collecting agencies in the Indus basin generally indicate to follow the WMO standards but we found inconsistencies in the use of precipitation measurement instruments and techniques. As the correction coefficients to account for wind-induced under-catch of precipitation depend on the type and orifice area of the precipitation gauge, incorrect gauge configuration information has consequences for the catch corrections. Although we tried our best to obtain the maximum possible information regarding the type and specs of precipitation gauges, we cannot exclude the chances of different precipitation gauges than the actual ones in some cases. However, we also think that the possibility of slight differences in gauge type will only have a small impact on the final results. The uncertainties resulting from spatial interpolation techniques described by Dahri et al. (2016) are equally applicable for this study as we followed their interpolation approach. Importantly, the cross-validation results infer high accuracy of the corrections and indicate excellent agreement between the adjusted precipitation and adjusted specific run-off at sub-basin scale.

| CONCLUSIONS
Reliable estimates of precipitation climatologies and amounts in the high-altitude Indus basin are seriously constrained by the quality and number of observed data (e.g., scarcity of in situ observations, measurement errors and space-time breaks). This study attempted to address these core issues by improved estimates of the precipitation Turc-Budyko representation of run-off ratio (Q/P) and aridity index (P/ET p ). The red triangles display estimates of unadjusted case or Dahri et al. (2016), black diamonds show estimates of Immerzeel et al. (2015) and blue circles indicate adjusted estimates under this study. The numbers refer to the sub-basins as given in Table 4 [Colour figure can be viewed at wileyonlinelibrary.com] measurement errors and integrating precipitation data from multiple sources with the net snow accumulations at major glacier zones. The study employed WMO recommended standard methods to adjust systematic errors in precipitation measurements. Simple methods to adjust net snow accumulation for the ablation losses and adjustment of river flows for the net mass balance contributions are introduced. Mean monthly adjusted and unadjusted precipitation observations and net snow accumulations are spatially interpolated using the KED interpolation scheme. Analysis of temperature variations with elevation and latitude revealed significantly different gradients for each month and substantial differences among the gradients at different locations for maximum and minimum temperatures. Hence, the use of a universal annual gradient or a time-independent gradient of mean temperature to estimate maximum and minimum temperatures or vice versa is a major source of uncertainty for the high-altitude Indus basin. The applied error-adjustments significantly increased the gauge-measured precipitation, which is in line with previous studies. The total bias between gauge-measured and erroradjusted precipitation ranged from 12 to 773 mm/year (2-182%) for various individual stations. The highest increments are computed for wind-induced under-catch of solid precipitation, particularly in higher-altitude areas and during winter months. The range of liquid precipitation under-catch is much smaller concentrating mainly in the low-altitude areas during summer monsoon. Similarly, notable increases varying from 0 to 1,000 mm/year (0-200%) are estimated for net snow accumulations. Precipitation increase at the basin (study area) scale is 21.3%, while at sub-basin scale it ranged from 6 to 77% with greater increments at higheraltitude areas and during winter months. Contrary to the general understanding, the contribution of net GMB to river flows is only marginal ranging from 0.5 to 6.1% of the observed flows. The highest contributions are revealed for the Chenab, Chitral, Shigar, Hunza, Astore and Gilgit basins.
The cross-validation results ( Figure 4) and the Turc-Budyko representation of the run-off ratios and aridity indices at sub-basin scale ( Figure 5) show that the adjusted precipitation amounts and distribution patterns derived in this study are more accurate than the unadjusted data and previous estimates. The catch corrections provided new insights in the magnitude and distribution patterns of precipitation implying potential hydrological implications for water resources assessment, planning and management. The actual precipitation is considerably greater than what has been previously thought. These increases are mainly realized in the higher-altitude areas of Chitral, Gilgit, Hunza, Shigar, Shyok and Astore basins. The study recognizes that the data quality-driven underestimated precipitation may be the major source of uncertainty in the water balance estimates in the high-altitude Indus basin. The improved climatologies of mean monthly precipitation developed in this study can be used for basin or sub-basinscale water balance studies and bias correction of gridded precipitation products, thereby paving the way for the development of an accurate, consistent and high-resolution gridded precipitation product for this highly under-explored region of the Indus basin.
Although our estimates of precipitation distribution can easily be regarded as much better than currently available estimates, the uncertainties elaborated at the end of the previous section recognize the need for further improvement. Further improvements can be achieved by calibration of the already installed precipitation gauges with the WMO recommended reference gauges and development of site and gauge-specific error adjustment models, use of observed data with better spatio-temporal coverage, use of daily or even sub-daily time steps, use of corresponding observed wind speed and temperature data sets, selection of any better spatial interpolation technique, accuracy assessment and precise determination of other components of the water balance to validate precipitation, and a better integration of precipitation data with mass balance data.