Correction and Informed Regionalization of Precipitation Data in a High Mountainous Region (Upper Indus Basin) and Its Effect on SWAT-Modelled Discharge

: The current study applied a new approach for the interpolation and regionalization of observed precipitation series to a smaller spatial scale (0.125 ◦ by 0.125 ◦ grid) across the Upper Indus Basin (UIB), with appropriate adjustments for the orographic effect and changes in glacier storage. The approach is evaluated and validated through reverse hydrology, and is guided by observed ﬂows and the available knowledge base. More speciﬁcally, the generated corrected precipitation data is validated by means of SWAT-modelled responses of the observed ﬂows to the different input precipitation series (original and corrected ones). The results show that the SWAT-simulated ﬂows using the corrected, regionalized precipitation series as input are much more in line with the observed ﬂows than those using the uncorrected observed precipitation input for which signiﬁcant underestimations are obtained.


Introduction
Water balance calculations and spatially distributed rainfall-runoff models require high resolution climatic datasets as a primary input. The quality of climatic forcing input is certainly the most important factor capable of influencing the simulation results [1][2][3][4][5][6][7], so that any errors in the input (climatic data) are amplified in the output (simulated hydrology) [8]. Among the climatic forcing data, precipitation and temperature are the most vital climatic variables as they directly influence the catchment discharge [9] through direct runoff, evapotranspiration losses or the snow and glacier melt contributions.
As precipitation is extremely variable in terms of spatial and temporal distribution, over mountainous catchments, it is extremely challenging to assess its true spatial distribution, especially when there is a limited spatial density of available gauging networks [10], and the then usually used "valley based" gauging networks are mostly unable to capture the orographic influences. The problem can be further exacerbated when the precipitation dataset has quality issues or temporal discontinuities. The hydrological investigations in such mountainous catchments may therefore find "water imbalances", with higher streamflow totals, in excess of precipitation-based estimates.
The precipitation data for the Upper Indus Basin (UIB) also suffer from these problems and the sparsely located, often lower-altitude, in-situ observational network, is unable to provide precipitation Being a high-mountain region, the UIB contains the greatest area of perennial glacial ice cover (~15,062 km 2 ) outside the polar regions, with 2174 km 3 of total ice reserves [25]. Some estimates [26] show even greater glacier cover, reaching up to 12% of the UIB (above 19,000 km 2 ). The altitude within the UIB ranges from as low as 455 m to a high of 8611 m and, as a result, the climate varies greatly within the basin [27].
The summer monsoon has little effect on the basin, as almost 90% is in the rain shadow of the Himalayan belt [22]. Except for the south-facing foothills, the intrusion of the Indian-ocean monsoon is limited by the mountains, so that its influence weakens northwestward [28]. Subsequently, the climatic controls in the UIB are quite different from that in the Himalayas on the eastern side. In fact, over the extent of the UIB, most of the annual precipitation originates in the west, resulting from the mid-latitude western disturbances, and mostly in solid form during winter and spring [19,20,29,30]. Occasional rains are brought by the monsoonal incursions to trans-Himalayan areas [20,29], but even during the summer months, the trans-Himalayan areas do not derive all precipitation from monsoon sources [31].
Climatic variables are usually strongly influenced by topographic altitude. Thus, the northern valley floors of the UIB are arid, with annual precipitation of only 100-200 mm. These totals increase to 600 mm at 4400 m elevation, and glaciological studies suggest annual accumulation rates of 1500 to 2000 mm at 5500 m altitude [20,32] (see Appendix A- Figure A1) The average snow cover area in the UIB changes from 10% to 70%, with a maximum of 70-80% in the winter snow accumulation period (December to February) and a minimum of 10-15% in the summer snow melt period (June to September) [27]. Stream flow is generated by the combination of the storm runoff in the lower part of the upper Indus basin and the snow and glacier runoff from the higher parts of the UIB [24,33]. Nearly half of the total annual flow in the Indus basin as a whole is contributed to by the UIB, with a 86-88% contribution during the summer season while only a 12-14% contribution during the winter season [34,35].
.   The stations operated by PMD have long term precipitation records (1947 to date), but were not corrected for under-catch and occasional gaps and temporal discontinuities in the records. The recent PMD-station records, although fairly consistent, are all from low altitude stations, unable to represent precipitation regimes of higher altitudes. The other 14 climate stations are fairly new and are installed specifically to cover higher altitude regions, but only provide data for a shorter time period from 1999-2008 (Table 1). For parts of the UIB outside Pakistan's boundary, including the Shyok basin and the UIB upstream of Kharmong (Figure 1), we used the APHRODITE-daily gridded precipitation dataset for Asia which is based on a dense network of rain gauges. [11]. This data has a spatial resolution of 25 km 2 and is available for the time period 1951-2007.

Observed Discharge Data
The daily river discharge and flow data in the study area are collected from the Water & Power Development Authority, Pakistan (WAPDA). Data for a total of 14 hydrometric stations is available in the UIB, out of which data for seven hydrometric stations (on the main river as well as its major tributaries) (Table 2, Figure 1) is used in this study.

Relevant Literature-Correction Methods
A range of methods (simple to very complex), have so far, been applied to improve the quality and coverage of precipitation data while trying to account as much as possible for the prevailing spatial and orographic variations. These methods can broadly come under two categories: an interpolation or regionalization of point rainfall measurements 2.
Both methodologies interpolation/regionalizing or correcting catchment-scale rainfall are extremely uncertain processes. For the former this is due to the spatial variability of rainfall fields and the complexities of orography, while for the latter (backward hydrology), uncertainties are faced due to the inherent nonlinearity in the streamflow-rainfall relationship.
To choose which of the two approaches to follow is difficult. Apparently there have been a lot more attempts to use interpolation/regionalization of point rainfall than the backward hydrology approach (e.g., References [36,37]).
For the first approach, the upward interpolation/regionalization of point rainfall measurements, many different techniques have been used in the past [38][39][40][41][42][43][44][45][46][47]. These techniques have been evaluated for performance, against data of different temporal scale and spatial coverage as well as for different regions, from having very simple and homogeneous terrain to those with highly complex and diverse topography. Nevertheless, their results and recommendations are as diverse as their application, which make a direct and conclusive comparison between the different variants of this method category difficult and impractical. In fact, it is very difficult to choose the method that could reproduce the climatic data distribution closest to reality in diverse catchment specific terrains [48], because the interpolation method which will best perform for one specific area, varies as a function of the area, terrain, the spatial scale desired for mapping [49], as well as the temporal duration and the nature of the climate variable to be interpolated [50].
Overall, the filling of spatial gaps through interpolation can be possibly done by three groups of techniques, i.e., empirical, statistical and geostatistical methods or function fitting [51].
Additionally, there are "other methods" which are specially developed for meteorological interpolation, using a combination of different methods, both deterministic and probabilistic [68]. An example of this type can be the "Meteorological Interpolation based on Surface Homogenized Data" (MISH), developed at the Hungarian Meteorological Service [69]. Though there is a range of methods in this category (simple to very complex and demanding) to choose from, the outputs from all of them carry huge uncertainties.
The second approach, i.e., "DHB" or Kirchner's Methodology, has mainly been used to infer the spatial and temporal patterns of evapotranspiration and precipitation at the gauged catchment scale, using measured streamflow fluctuations as input. This method has been reported to be effective in estimating catchment-averaged precipitation (e.g., References [21,36,37,70,71]), but as it may not give distributed precipitation, any explicit regionalization of precipitation fields may need extra efforts or modeling [72].
Any of these techniques or strategy, which can be used for correction and regionalization of data, need calibration and validation by means of historical information, [49], directly or indirectly, by evaluating the results and outputs of spatially distributed hydrological models.
In the case of the interpolation method, when the study area has a dense enough coverage of available observations, the data is divided in to two sets, with one set used for interpolation and the other one for independent validation. When this is not the case and there are considerable spatial gaps in the data points, the comparison is usually done through cross-validation [73]. In cross-validation, data at a gauge-point is removed temporarily, one at a time, and re-estimated from the remaining data. The estimated values are checked against the observed values to evaluate the accuracy of the interpolation methods. These validation techniques can be further supplemented by an indirect validation, wherefore streamflow observations are used as a reference for evaluating the results and outputs of hydrological simulation models. This kind of validation needs more efforts, because the correction/regionalization method needs to be combined with a spatially distributed hydrological model to evaluate the results of both combined [71,[74][75][76][77]. Nevertheless, in some cases (such as for the method proposed in this study) this combined approach becomes the only viable validation option. Although this method is more complex in terms of efficiency assessment, it can provide more options to assess the results over a wider spatial and altitudinal range, depending on the availability of stream flow data [21].

Method Used in the Present Study
The informed spatial regionalization of the precipitation data in the UIB involves basically three main processes: (1) Pre-processing of the raw precipitation data (quality check & correction of systematic errors), (2) estimation of a catchment-specific orographic correction factor (OCF) or precipitation laps rate (PLR) by "Doing Hydrology backward", and (3) stepwise interpolation/regionalization (OCF adjustment of observed precipitation at catchment mean elevation, followed by simple kriging and finally an OCF readjustment of the interpolated data to target grid-average/point elevation). These three steps are described in detail in the subsequent sub-sections.

Correction of Systematic Errors
The first step consists of the correction of systematic errors in the raw uncorrected precipitation data available for the UIB. This task serves to remove any systematic errors in precipitation measurements and to correct the possible precipitation under catch during snowfall, particularly, under windy conditions. For this purpose, different methods recommended by the World Meteorological Organization (WMO) were reviewed, and the method of Reference [78] was selected as it requires less observed parameters and has also previously been applied in the study region. It should be noted though that, prior to application of the corrections for systematic errors, the observed precipitation dataset was checked for inhomogeneity and outliers. This correction method employs Equations (1) Water 2018, 10, 1557 7 of 23 and (2) to account for wind-induced errors, wetting losses, evaporation losses and trace amounts. The equations suggested by the method in References [78,79] are as follows: where Pc is the 'corrected' precipitation; Pm, the measured gauge value; ∆Pw, the wetting losses; ∆Pe, evaporation losses; ∆Pt, trace amount; and K, the adjustment coefficient due to wind-induced errors, for which CR is the catch ratio (%), defined as a function of wind speed. The values of ∆Pw, ∆Pe, ∆Pt and CR suggested by References [78,79], used in this study, are given in Table 3. For allocating values for ∆Pw, ∆Pe, ∆Pt and CR in the calculations, the precipitation was considered during (1) the winter months (DJF) as Snow, (2) the summer months (JJO) as Rain, and (3) spring and autumn (MAM and SON) as Mixed.

Estimating the "Orographic Correction Factor" (OCF) by Doing Hydrology Backward (DHB) plus Interpolation (Informed Regionalization (IR))
• General approach The "Orographic correction factor" (OCF) is calculated based on the rearranged "hydrological/water balance equation" (Equations (6) and (7)). The calculation utilizes five data variables, out of which data for three variables: mean annual catchment precipitation, mean annual catchment discharge, and catchment and gauge point elevation were available, whereas the other two variables required, catchment mean annual change in glacier storage (mm) and catchment mean annual actual-evapotranspiration, were estimated based on the relevant literature and gridded data products. The final OCF is through an adjusted variant of the "equation-based version", which was computed based on the hydrological modeling calibration, until the simulated mass balance came to an acceptable match to the observed mass balances.

•
The hydrological/water balance equation The hydrological/water balance equation [80] is as follows: where Q t is the total volume of water discharging from a catchment (per specified time period), and is equal to the volume of water entering the catchment as true precipitation (P) minus a change in storage and losses from the system. The latter comprise of losses by actual evapotranspiration (ET), aquifer/ground-water recharge (Gw) and losses or gains of glacier ice volume (∆g).
As the losses Gw, may be minimal in comparison to the total discharge especially for large mountainous catchments and at an inter-annual time step, the above equation can be simplified to: Water 2018, 10, 1557 By assuming the true aerial precipitation P true to be equal to the observed precipitation P obs plus the orographic under/overestimation OCF, one gets: Basically, the OCF T represents the under-or over estimation by the low-altitude average observed precipitation of the true areal precipitation over the gauged catchment.
The OCF per unit elevation, OCF plapse , is found by dividing OCF by the difference ∆h in mean elevation of the catchment and the observation network, i.e., by arranging Equations (6), one gets One can also define the orographic correction multiplicative factor, OCF multiplicative , by dividing the true precipitation P true by the observed precipitation P obs : In the presence of the water discharge data from the system (Q), observed precipitation (P obs ) and elevation (∆h), only estimates of glacier mass balance (∆g) and evapotranspiration (ET) are required to calculate the OCF for each gauged catchment. The values for these parameters are estimated based on available literature values for gridded data sets for glacier mass balance and actual evapotranspiration, as detailed in the subsequent sub-chapters.

•
Glacier mass balance (∆g) estimates Estimating the glacier storage change (∆g) in the different catchments of the UIB is a very difficult and uncertain task, because no consensus could be found within the available literature on either the amount of decrease in most of the southern and eastern parts or the possible increase (or decrease) in the northern or western parts of the HKH region.
Much of the UIB spans over different mountain ranges, out of which the watersheds of Hunza, Shigar, and Shyok span the Karakoram mountains; Gilgit watershed covers the north-eastern Hindu Kush, while four watersheds (Astore, Shingo, Zanskar, and remaining part of the UIB upstream of Kharmong) are located across the Western Greater Himalayas. Therefore, the mass balance of the glaciers in these mountain ranges and the watersheds spanning over them is as diverse as their hydro-climatic regimes are different and hard to find out.
Although the decreases or increases suggested by different studies may differ notably, many of the recent studies are consistent in pointing out a decreasing trend in the eastern and central Himalayas and a steady or slightly increasing trend in the Karakorum and other adjoining ranges in the north and west. These claims are supported by other hydro-climatic factors, such as an increase in observed precipitation for winter and summer [81,82], the Karakorum anomaly validated by Reference [83], with the findings of a positive mass balance for the central Karakoram glaciers and less summer flows, in spite of a precipitation increase [82].
The snow cover change over different catchments in the UIB (Figure 2), derived based on statistics available at "HKH Snow Cover-web application" hosted by ICIMOD, is also in conformity with the claims of a stable to positive glacier mass balance in the Karakorum. Other recent studies [32,[82][83][84][85][86][87][88] indicate similar trends of glacier mass balance in the region. These trends of glacier mass balance are summarized for a few recent studies in Table 4.  Despite these agreements on the trends of net glacier mass change, there are huge differences between the studies regarding their absolute amounts. These uncertainties are further exacerbated by the fact that different authors used different source data, different assessment periods and different procedures. To avoid some of these complexities, one of the recent studies, i.e., Reference [83], which is also in line with others (e.g., References [82,87,88] etc.), is selected as a reference for glacier mass balance in different parts of the UIB, with the major results as outlined below: (a) For the western and northern parts of the UIB, namely Hunza and Shigar, the glacier mass balances are assumed to have a net positive trend [82,88] and are allocated a value of +0.09 mwe year −1 as proposed by Reference [83] for the western Karakorum. (b) The Gilgit watershed, though spanning over the northern end of the Hindukush range (bordering the western Karakoram) is assigned the same value of +0.09 mwe year −1 , because, as reported by Reference [82], this watershed has more similarities with the climatic regime of the Karakoram region (Hunza basin), rather than the rest of the Hindukush area. (c) For the Shyok river basin, a snow cover change of +0.11 mwe year −1 is assumed based on the value given by Reference [83] for the east Karakorum.   Despite these agreements on the trends of net glacier mass change, there are huge differences between the studies regarding their absolute amounts. These uncertainties are further exacerbated by the fact that different authors used different source data, different assessment periods and different procedures. To avoid some of these complexities, one of the recent studies, i.e., Reference [83], which is also in line with others (e.g., References [82,87,88] etc.), is selected as a reference for glacier mass balance in different parts of the UIB, with the major results as outlined below: (a) For the western and northern parts of the UIB, namely Hunza and Shigar, the glacier mass balances are assumed to have a net positive trend [82,88] and are allocated a value of +0.09 mwe year −1 as proposed by Reference [83] for the western Karakorum. (b) The Gilgit watershed, though spanning over the northern end of the Hindukush range (bordering the western Karakoram) is assigned the same value of +0.09 mwe year −1 , because, as reported by Reference [82], this watershed has more similarities with the climatic regime of the Karakoram region (Hunza basin), rather than the rest of the Hindukush area. (c) For the Shyok river basin, a snow cover change of +0.11 mwe year −1 is assumed based on the value given by Reference [83] for the east Karakorum. (d) For the region covered by the Astor basin and the parts of UIB downstream of Kharmong, except the aforementioned tributary catchments, the selection is not straight forward, as some of the above studies claim a constant/slight increase in snow cover area in this region (e.g., Reference [82]), while others suggest a negative mass balance. Therefore, for the current study, the glacier mass balance in this region is taken as neutral, with no increase or decrease. (e) For the parts of UIB, upstream of Kharmong, which have experienced a consistent decrease in snow/glacier cover, according to most of the literature, a value of −0.45 mwe year −1 is assigned, based on the mean mass balance value [83] for the Spiti-Lahaul region.
The glacier mass balance amounts selected above for the different regions are converted from meters water equivalent per year (mwe year −1 ) to glacier storage changes in millimeter depth (mm) ( Table 5), before using them in the calculation for the orographic correction factor.
• Evapotranspiration (ET) estimates: Observed actual evapotranspiration (ET) data are a very rare commodity, even in easy-to-access areas. This holds also for the UIB, where there are not much direct measurements available [72]. For ET-estimations, especially in the UIB and elsewhere too, most of the previous studies relied on either temperature-based empirical relationships or, in more recent times, on calculated algorithms utilizing data with different combinations of various ground-or sensors-based observed hydro-climatic variables [90][91][92][93][94]. A few studies suggested to address the uncertainties in the different ET-estimates by blending them to acquire an average trend (e.g., Reference [72]).
Despite some improvements in the estimation techniques or the data used in these studies, the uncertainties in ET-estimates are huge, especially for areas like the UIB, where the other available data also have huge uncertainties. For the UIB, most of the gridded evapo-transpiration products show very high average annual ET values, which appear unrealistic when compared to the observed annual precipitation over the UIB of 350-400 mm year −1 (~367 mm year −1 based on the observed and APHRODITE precipitation data, given in and Table 7 later). The mean annual ET values, proposed by most of the available gridded products such as the "Global Average Annual Evapotranspiration, 1950-2000", data hosted at Databasin.org [95], have very high average ETs over different catchments of UIB (~300-400 mm year −1 ). Similarly, the ET values proposed by Immerzeel et al. [72], based on four different available ET products, have a mean annual value of 359 ± 107. These proposed mean annual ET values are probably too high to leave any water left for runoff generation.
The reason is that in the high-elevation mountainous catchments of the HKH, evapotranspiration only plays a minor role in the overall water balance and is generally less than 10% of the total hydrological budget [96]. This also applies to the UIB where only the lower-altitude area may have higher ET contributions, while the high elevation sub-catchments may have only minor annual ET due to their low average temperatures. For example, in the Hunza basin of the UIB, Garee et al. [97], reported a model ET equal to 18% of the precipitation, indicating even lower ET in the higher altitude catchments. Therefore, in UIB, which is also a high mountainous basin and has a mean annual water discharge of around 462 mm year −1 , the ET in UIB should not amount to such high values, almost equal to the observed runoff.
To solve these issues, ET products with relatively lower values were evaluated. One such ET data is the Esri_hydro "Average Annual Evapotranspiration" [98], which is based on the MOD16 Global Evapotranspiration Product, and is derived from MODIS-satellite imagery by a team of researchers at the University of Montana, has comparatively lower ET-estimates, which better match the recommended values of Bookhagen and Burbank [96] or the ET-values reported by Garee et al. [97]. The MOD16-ET-data have a good resolution of 1 km 2 and are available over the period 2000-2011. For this reason, they are used here as reference-ET in the OCF-calculation, and their values are also listed for the different UIB-catchments in Table 5.

Regionalization Procedure-Step-Wise Interpolation
The gauge-station precipitation records and the APHRODITE-precipitation were processed separately. The gauge-station-observed precipitation time series, after correction for systematic errors and quality check, was interpolated to 0.125 • by 0.125 • grid in three steps.
In the first step of this "regionalization/stepwise interpolation" process, all the point-observed data are adjusted for the mean catchment elevation as: where P target is the precipitation at target elevation (mm); P gauge is the precipitation recorded at the gauge station (mm); EL target is the elevation at the target point/grid; EL gauge is the elevation at the gauge station; wd is the average annual number of wet days; and OCF plapse (Equation (8)) is the orographic correction factor, in terms of the precipitation lapse rate, for the catchment.
In the second step, the data adjusted at the mean catchment elevation is then interpolated using "Simple Kriging" [59], to a 0.125 • by 0.125 • grid as well as SWAT-Model sub-basin centroids.
In the third step, the interpolated data are readjusted from the mean catchment elevation to the grid elevation, using Equation (9).
For the APHRODITE data, the correction for the orographic effect was done using the multiplicative correction factor OCF multiplicative (Equation (8)), prior to interpolation using "Simple Kriging" [59], to a 0.125 • by 0.125 • grid as well as SWAT-Model sub-basin centroids, as discussed below.

• SWAT Model Description
The Soil and Water Assessment Tool (SWAT) is a hydrological model developed for the US Department of Agriculture (USDA), Agricultural Research Service (ARS) by Dr. Jeff Arnold and collaborators. The SWAT-model is a continuous time (long-term yield), process-based semi-distributed model, capable of simulating hydrological processes in river basins/watersheds, based on specific information pertaining to the watershed, such as weather/climate, topography, soil properties, land cover, land use and management practices [99,100]. In the SWAT-model, a main river basin or watershed is partitioned into several sub-units called sub-basins draining the tributaries into the stream network and the river-system. These sub-basins are further divided into a series of smaller units, the co-called hydrological response units (HRUs), which are spatial uniform units, each representing a unique combination of soil, land-use and slope. The calculation and simulation of the various hydrological components is based on the solution of the fundamental water balance equation (Equation (3)), wherefore these components which may include, in addition, sediment yield and agricultural nutrients, are first evaluated for each HRU and then routed and aggregated for the subbasin and finally for the watershed.

• Model Calibration and Validation setup
The SWAT model was calibrated and validated against daily discharge individually for each of its five major tributaries (Hunza, Gilgit, Astor, Shigar and Shyok rivers), for parts of UIB (except the tributaries) inside Pakistan's boundary and for the UIB (situated in India China and Nepal) covering the area upstream of the Kharmang gauge station. In cases of catchments, where inflow from the upstream catchment had to be accounted for, the observed discharge was used as inflow.
The Sequential Uncertainty Fitting SUFI-2 algorithm [104] of the SWAT-CUP program [104] was used for parameter optimization during the calibration process.
For the performance evaluation of calibration/validation results, the "goodness of fit" statistics including the coefficient of determination (R 2 ), Nash Sutcliff efficiency (NS) and Percentage bias (PBIAS) are used, which assess the simulated hydrological responses against the observed flow data.

Construction of Orographically-Corrected Precipitation Datasets
The selected and finalized values for the glacier mass balances ∆g and the actual evapotranspiration ET (Table 5) are utilized in Equations (7) and (8) to derive the two kinds of catchment-specific orographic correction factors (OCF), namely, the multiplicative correction factor OCF multiplicative and the additive correction factor OCF lapse (representing the precipitation lapse rate per 1000 m elevation change). The former are derived for the uncorrected raw precipitation for both gauge-and APHRODITE-data sets, while the latter is calculated only for the gauge station records, after correction for systematic errors.
The additive correction factor OCF lapse is applied at all the elevation ranges, wherefore it is assumed that the precipitation increases uniformly with elevation and the correction factor is constant throughout the year. To further improve the regionalization, and with more data available, it would also be possible to apply a range of catchment-or season-specific OCF's to generate the desired vertical or temporal regimes of precipitation to match any empirical distributions of precipitation intensities in different elevation zones or seasons.

True Areal Precipitation and OCF's
The details of the annual input data and the results obtained for the different OCF's for the various UIB catchments are listed in Tables 6 and 7. The findings are in total conformity with conclusions of many previous studies, which claim that the gauge-based data as well as the remotely sensed precipitation products are unable to represent the true areal precipitation in the UIB [11][12][13][14].  Notes: a applicable to raw gauge precipitation records, b applicable to gauge precipitation records already corrected for systematic errors, c applicable to gauge precipitation records already corrected for systematic errors as additive factor per 1000 m, * UIB upstream of Kharmong including Shingo, Zanskar; ** UIB downstream of Kharmong without main tributary catchment, + APHRODITE data is used only for Shyok and Kharmong basins, ++ Observed gauge station records, corrected for systematic errors through Equation (1), +++ True areal precipitation estimated based on Equation (5).
• Spatial distribution of orographically-corrected precipitation datasets The final gridded precipitation generated after the step-wise interpolation to the 0.125 • × 0.125 • grid by simple Kriging, followed by the appropriate elevation correction with OCF lapse at the target points is shown in the map of Figure 3a. Obviously the application of the elevation OCF lapse , while remaining true to the originally estimated mean precipitation for each catchment (Table 6), has added some diversity, as it produces a spatially distributed precipitation over UIB and its catchments with elevation-dependent lapse gradients. For comparison, the observed/APHRODITE interpolated precipitation distribution as well the difference between the two are shown in Figure 3b,c, respectively, with the latter witnessing clearly the underestimations of the observed/APHRODITE precipitation in most areas of the UIB.
Overall, the corrected and regionalized precipitation in Figure 3a, reveals outstanding patterns of spatial and orographic distributions. Our results match well with those of recent available studies of precipitation with respect to intensities, horizontal and vertical distribution as well as regional trends in the UIB [23], although our precipitation values are slightly lower as [23], because we have assumed a smaller ET rate over the UIB. The spatial distribution and trends obtained here are also in conformity with others (e.g., References [19,20,24,30]), such that over the extent of the UIB, the highest average annual precipitation is found in the west, due to the prevailing midlatitude western disturbances. The monsoonal contributions over the UIB are also in accordance with results of (e.g., References [20,24,28]), such that they act mostly in the southern fringes of the western Himalayas in the UIB, but are waning in the north and west (Figure 3a). the UIB, but are waning in the north and west (Figure 3a). Table 6 shows that the mean OCF-corrected precipitation over the UIB for the 1999-2008 period has a value of ~608 mm/year, which is considerably higher than the average uncorrected gauge station and APHRODITE measured annual precipitation of 367 mm/year, (i.e., the latter is underestimated by ~166%). For some catchments in the UIB (Gilgit and Shyok), these underestimations of the true precipitation are even in excess of 300%. The strong relation between the topographic altitude and precipitation amount [20], is also well accounted for by our results, with very high precipitation over the elevated zones of the Gilgit and Astore basins, and comparatively weaker orographic influences in the rain-shadow regions in the north or east. Resultantly, the north-western parts of the UIB, which lie inside Pakistan's boundary, show the highest mean annual precipitation, especially over the Gilgit and Astor river catchments, of above 2000 mm/year. There is a gradual decrease in the mean annual precipitation to the north of these two catchments, with Hunza and Shigar having lower precipitation. A similar decreasing  Table 6 shows that the mean OCF-corrected precipitation over the UIB for the 1999-2008 period has a value of~608 mm/year, which is considerably higher than the average uncorrected gauge station and APHRODITE measured annual precipitation of 367 mm/year, (i.e., the latter is underestimated by~166%). For some catchments in the UIB (Gilgit and Shyok), these underestimations of the true precipitation are even in excess of 300%.
The strong relation between the topographic altitude and precipitation amount [20], is also well accounted for by our results, with very high precipitation over the elevated zones of the Gilgit and Astore basins, and comparatively weaker orographic influences in the rain-shadow regions in the north or east. Resultantly, the north-western parts of the UIB, which lie inside Pakistan's boundary, show the highest mean annual precipitation, especially over the Gilgit and Astor river catchments, of above 2000 mm/year. There is a gradual decrease in the mean annual precipitation to the north of these two catchments, with Hunza and Shigar having lower precipitation. A similar decreasing precipitation gradient is witnessed in the west-east direction of the UIB, so that the lowest mean precipitation is witnessed in the easternmost parts of the basin.

Validation of the Corrected Precipitation against SWAT-Simulated Discharge
To validate the final product of the proposed precipitation correction and regionalization method, the SWAT hydrological model was forced with the two different precipitation data series: (1) observed, gauge-and (2) orographically-corrected & regionalized precipitation, both for the time period 1999-2008.
The goodness of fit statistics of the SWAT-simulated discharges for the two categories of precipitation input data are listed in Table 8, while the comparisons of simulated and observed discharges (hydrographs) are presented in Figure 4. That is why, when forced with corrected and regionalized precipitation, the SWAT-simulated daily discharges match the observed discharges much better, and this holds for all the monitoring points across the various tributary catchments (Figure 4).
The goodness of fit statistics (Table 8) also show drastic improvements over those of the simulations forced with uncorrected gauge station precipitation series. For all the monitoring points in the UIB, the R 2 are above 0.69 and go as high as 0.89. The NS-values are also on the higher side and range between 0.69 and 0.89. Overall, a R 2 of 0.86 and a NS of 0.85 is obtained at the basin outlet, which can be considered as very good results, especially for hydrological simulations carried out at a daily time step. Furthermore, the PBIAS-values of the simulated discharges are also much lower than those of the first SWAT-simulations category, with maximum positive and negative biases of 12.4% and −19.7%., respectively.
These results verify clearly the validity of the adopted methodology of precipitation correction and informed regionalization by accounting for orographic influences and regional patterns at the gauged catchments scale for use in hydrological modeling studies.

Merits and Limitation of the Proposed Method and Regionalized Data
We will assess our precipitation data in comparison of the available precipitation products as well as check the comparative strengths and weaknesses of our proposed method.
Currently, various gridded precipitation data sets are available for HKH. They include satellite data products, such as the TRMM 3B42 product, rain-gauge based collections (i.e., APRHODITE, CRU, GPCC), combined satellite and rain gauge climatology (GPCP), reanalysis products (i.e., ERA-Interim) etc. Although all these data sets are reasonably suitable for hydro-climatic studies at a regional scale; owing to the complex orography of the UIB region and to the co-action of different hydro-climatic regimes (which affect the amounts, spatial patterns and the seasonality of precipitation), neither the sparse observed station data (or the gridded data products based on them) nor the sensors-based data, fully represents the precipitation regime of the region [11][12][13][14].

Merits and Limitation of the Proposed Method and Regionalized Data
We will assess our precipitation data in comparison of the available precipitation products as well as check the comparative strengths and weaknesses of our proposed method.
Currently, various gridded precipitation data sets are available for HKH. They include satellite data products, such as the TRMM 3B42 product, rain-gauge based collections (i.e., APRHODITE, CRU, GPCC), combined satellite and rain gauge climatology (GPCP), reanalysis products (i.e., ERA-Interim) etc. Although all these data sets are reasonably suitable for hydro-climatic studies at a regional scale; owing to the complex orography of the UIB region and to the co-action of different hydro-climatic regimes (which affect the amounts, spatial patterns and the seasonality of precipitation), neither the sparse observed station data (or the gridded data products based on them) nor the sensors-based data, fully represents the precipitation regime of the region [11][12][13][14]. This has led many researchers to find ways to assess methods for precipitation correction that may lead to a more realistic water balance [21] in many basins and have also compelled a number of hydrological studies in the UIB region to use, in addition to the observed station data, a variety of other reference climate data from different sources, either directly or with prior modifications and adjustments (e.g., TRMM Data [22]; modified APHRODIT [23], modified WFDEI data [13] or modified observed data through reverse hydrology [72] etc.). All these data sets have been produced through reasonably complex methodology or possibly have under or overestimated precipitation magnitudes. For example, the precipitation estimates by Reference [72] have utilized evapotranspiration values which appear to be higher than the realistic (or reported) limits.
To avoid these downsides, we have tried to devise a method which may not only be less complex, but also produce precipitation data capable of capturing the prevalent spatial patterns and orographic influences in UIB, at least at gauged catchment scale. The proposed method is not only guided by realistic or reported glacier mass-balance, evapotranspiration and orographic influences, but also producing precipitation estimates that are able to close the water balance in UIB.
Although our results clearly suggest the validity of the proposed method, the method is not without limitations. The spatial patterns of precipitation are based on a fixed correction factor for each gauged catchment, and hence, are unable to provide explicitly distributed precipitation, as the single orographic correction factor applied per catchment cannot compensate for the variation in orographic influences inside a catchment. Similarly, the estimated precipitation inherited uncertainties attached to the glacier-mass balance estimates as well as ET data used.

Conclusions
In this study, we addressed the underestimation and low representation of the high-altitude precipitation by the available gauge-based records by doing hydrology backwards in conjunction with available cryosphere and hydro-climatic information. Despite certain uncertainties, our precipitation estimates are not only accounting for the orographic influences, but also for the glacial mass balance across different catchments of the UIB. The corrected and regionalized data is therefore making it possible for hydrological investigations to properly close the water balance, which is unlikely to be achieved with the available observed gauge-based or gridded precipitation datasets.
The estimations of the horizontal and vertical distributions as well as the magnitudes and intensities of the precipitation achieved by our correction and informed regionalization approach are matching well with those reported in the literature. Furthermore, the SWAT-hydrological simulations (doing hydrology backwards) using the corrected precipitation data match the observed flow in the UIB even better than our expectations.
Thus, whereas the SWAT-hydrological simulations using uncorrected observed precipitation show poor fit with the observed discharges, with values for R 2 and NSE not greater than 0.77 and 0.41, respectively, and huge underestimations of the flows across all the investigated catchments of the UIB, with PBIAS-values of up to −72.87%, the situation is much improved when using the corrected precipitation datasets in the model. In these cases, simulated flows match the observed flows very well, with R 2 ranging between 0.69 and 0.89, with a similar range for the NSE and, last but not least, small PBIAS-values ranging between 12.4% and −19.7%.
The results of the present study show that major improvements in rainfall estimations in poorly gauged, high mountain regions, like the UIB, can be achieved by combining classical orographic correction methods with knowledge of the regional hydro-meteorology and glacier mass balance at the gauged catchment levels (in case of glaciated catchments) and validating such a methodology by an additional hydrological model, in order to remove inconsistencies in and to close the hydrological water balance (doing backward hydrology). This multistep approach is not only less demanding in terms of computational or human resource requirements than the methods involving advanced atmospheric physics or geostatistics, but can considerably improve the quality and representativeness of the precipitation data at a gauged catchment scale.