Study of the Spatiotemporal Characteristics of Meltwater Contribution to the Total Runoff in the Upper Changjiang River Basin

Melt runoff (MR) contributes significantly to the total runoff in many river basins. Knowledge of the meltwater contribution (MCR, defined as the ratio of MR to the total runoff) to the total runoff benefits water resource management and flood control. A process-based land surface model, Noah-MP, was used to investigate the spatiotemporal characteristics of MR and MCR in the Upper Changjiang River (as known as Yangtze River) Basin (UCRB) located in southwestern China. The model was first calibrated and validated using snow cover fraction (SCF), runoff, and evapotranspiration (ET) data. The calibrated model was then used to perform two numerical experiments from 1981 to 2010: control experiment that considers MR and an alternative experiment that MR is removed. The difference between two experiments was used to quantify MR and MCR. The results show that in the entire UCRB, MCR was approximately 2.0% during the study period; however, MCR exhibited notable spatiotemporal variability. Four sub-regions over the Qinghai-Tibet Plateau (QTP) showed significant annual MCR ranging from 3.9% to 6.0%, while two sub-regions in the low plain regions showed negligible annual MCR. The spatial distribution of MCR was generally consistent with the distribution of glaciers and elevation distribution. Mann-Kendall (M-K) tests of the long-term annual MCR indicated that the four sub-regions in QTP exhibited increasing trends ranging from 0.01%/year to 0.21%/year during the study period but only one displayed statistically significant trend. No trends were found for the peak time (PT) of MR and MCR, in contrast, advancing trend were observed for the center time (CT) of MR, ranging from 0.01 months/year to 0.02 months/year. These trends are related to the changes of air temperature and precipitation in the study area.


Introduction
Precipitation in the form of rainfall can trigger runoff generation process when it reaches the land surface and a specific threshold is satisfied, e.g., the rainfall intensity exceeds the infiltration capacity of the soil or rainfall occurs in an area with saturated soil.In contrast, precipitation in the form of snowfall is stored on the land surface and does not release water until atmospheric conditions reach a specific combined threshold and melting occurs.Melt runoff (MR), i.e., the runoff that is generated by melt of seasonal snow packs and glaciers, contributes significantly to overall runoff in many basins worldwide, especially during warm, dry summers.Runoff constitutes the primary freshwater resources in many areas [1,2], such as the Indus Basin [3], the Tien Shan Basin [4], the Colorado River Basin [1], and others.However, MR can also result in flood events during the spring and summer [5,6].These events complicate reservoir operation and flood control.The response of snow/glaciers to climate change (mainly changing temperature and precipitation) which affect the MR and the entire hydrological cycle is becoming more complicated [7,8].
The contribution of MR (MCR, expressed as the ratio of MR to the total runoff) to the total runoff is a commonly adopted indicator used to represent the importance of MR to total runoff [3].Although considerable challenges exist due to the complexity of hydrological processes, researchers have studied these processes in depth and proposed several analysis methods to investigate MCR at large scales [9].These approaches can be divided into two categories: observation-based methods and model-based methods.Observation-based methods use the ground measured or remotely sensed data to obtain runoff mass [10,11] or chemistry [12,13] information to quantify MCR based on certain assumptions [9].Observation-based methods are straightforward; however, the availability, accuracy, and representativeness of data are not always guaranteed, especially in remote regions or over large areas.By comparison, model-based methods use numerical models ranging from simple empirical models to sophisticated process-based models that simulate the key variables related to MCR.Model-based methods are widely used to study MCR because they can reproduce major hydrological components after proper calibration and validation, and they provide detailed information that is not available by observations [14][15][16].Another obvious advantage of the model-based approach is the ability to project future MR scenarios under changing environmental conditions [3,4,17].Among these models, Land Surface Models (LSMs) have been developed to understand the water, energy and carbon cycles of Earth [18].Many models emphasize on snow/glacier processes due to their important roles in weather and climate [19].Therefore, snow-related processes including accumulation, sublimation and melt are generally represented using physically based schemes in LSMs by explicitly considering various controlling factors, e.g., climate, topography, vegetation, etc.Therefore, LSMs are preferred for investigating MR and MCR, especially over large areas [1,15,20].
The Changjiang River (A.K.A. Yangtze River), originating in the glaciers of the Qinghai-Tibet Plateau (QTP) and flowing easts to the East China Sea, is the longest river in China.The river supports numerous anthropogenic activities in China; therefore, flood control and water resources management have long been major concerns of local authorities.Researchers have noted the importance of meltwater to the Changjiang River and performed related studies (e.g., [20][21][22][23][24]).However, these studies are limited to the source region (refers to the area above the Zhimenda runoff gauging station) of the Changjiang River Basin.Because the Changjiang River Basin covers different topographic and climate zones, the meltwater contributions to runoff over the large area of the Changjiang River Basin are unknown.Climate change is expected to affect the timing and magnitude of runoff of the Changjiang River, and meltwater is an important factor in these processes [25].Moreover, in addition to the Three Gorges Reservoir, many reservoirs are being constructed or planned on the main stream and tributaries of the Changjiang River.Understanding the spatiotemporal characteristics of MR can help improve water resources management, flood control and climate change adaptation in the Changjiang River Basin.
The overall objectives of this paper are to investigate the spatiotemporal characteristics of MR and associated MCR in the Changjiang River Basin using the process-based land surface model, Noah-MP.We restrict our study to the upper reach of the Changjiang River Basin (UCRB) because it is mainly located in the QTP region where snow, glaciers and meltwater tend to have larger effects on runoff, and the area is more sensitive to climate change [26].Additionally, reservoirs being planned and constructed are mainly located in this area.

Study Area
UCRB covers a large area of approximately 1 million km 2 from the source region on the QTP to Yichang in Hubei Province.The total length of the main stream is 4504 km (Figure 1).The distribution of glaciers extracted from "The second Chinese glacier inventory" [27] is also shown in Figure 1.The corresponding glacierized area is approximately 1687.2 km 2 .Due to the large geographic extent of UCRB (24.E-111.5 • E), there are large differences with respect to topographic and hydrometeorological factors.For example, the averaged elevations of the source region and Jialing river basin (a subregion of the UCRB, see below for details) are 4760 m and 1265 m, respectively.By comparison, the average annual air temperatures in these two regions are −5 • C and 16 • C, respectively.Such differences lead to diverse hydrometeorological conditions in the UCRB.To better analyze the various hydrometeorological processes, we divided the basin into 10 sub-regions based on the drainage network and runoff gauging stations.Note that the sub-regions defined in this paper are slightly different from those listed in the hydrological yearbook released by the Bureau of Hydrology, Ministry of Water Resources of China.Table 1 lists the characteristics of each sub-region, including the name, area, land cover and runoff gauging stations.Seven atmospheric variables required by the model were retrieved from the China Meteorological Forcing Dataset (CMFD [29]) (precipitation, downward short wave radiation, downward long wave radiation, wind speed at 2 m, air temperature at 2 m, air pressure and specific humidity at 2 m).The dataset was produced by merging a variety of data sources including ground observations from the China Meteorological Administration (CMA), Princeton reanalysis data, Global Land Data Assimilation System (GLDAS) dataset, etc. CMFD data were available between 1979 and 2015 with spatial and temporal resolutions of 0.1 • and 3 h, respectively.

Land Data
Land cover and soil types are two important land surface characteristics that determine model parameters.The land cover dataset is the 0.5 km MODIS-based Global Land Cover Climatology dataset released by Broxton et al. [28] based on 10 years (2001-2010) of MODIS land cover (MCD12Q1, Collection 5.1) products.As reported by Broxton et al. [28], it provides more realistic land cover classification over long-term period.The soil dataset we chose is "the Soil Database of China for Land Surface Modeling" developed by Shangguan et al. [30].It was derived from 8979 soil profiles and the Soil Map of China (1:1,000,000), which provides soil types at a spatial resolution of 30 .

Observed and Remote Sensing-Based Data
We used the observed data including runoff, snow cover fraction (SCF), and evapotranspiration (ET) data to calibrate and validate the model.

Runoff
Runoff data are collected and processed by the Bureau of Hydrology of the Changjiang Water Resources Commission.Based on the data availability, we chose six runoff gauging station in this study (Figure 1).Among these stations, the Yichang station is located at the outlet of the entire UCRB; however, runoff at the station is strongly affected by the reservoirs including the Three Georges Reservoir.To eliminate this effect, we used the computed natural runoff series at the Yichang runoff gauging station in this study.The natural runoff is computed by taking into consideration the reservoirs and anthropic activities.The remaining five runoff gauging stations control the corresponding sub-regions (Table 1).We adopted the measured runoff directly at these stations for calibration and validation, since no computed natural runoff available.
The temporal resolution of the original data varies from several hours to several days.Therefore, we derived monthly runoff using the time-weighted average method.

Snow Cover Fraction
We adopted the MODIS/Terra monthly SCF product (MOD10CM Version 6) to evaluate the model performance in modeling snow cover dynamics.The 0.05 • estimation of the SCF is generated by spatially and temporally aggregating the daily normalized difference snow index (NDSI) at a 500 m spatial resolution [31].MOD10CM provides information from March 2000 to present; however, there is no data available for four months (August 2000, June 2001, March 2002 and December 2003).

Evapotranspiration
Zhang et al. [32] derived an ET dataset (ET_Zhang) based on a remote sensing ET algorithm, providing global monthly ET records from 1983 to 2006 at an 8 km spatial resolution.Studies using the ET_Zhang dataset in the Changjiang River Basin have confirmed its accuracy [32,33].In this paper, we used the ET_Zhang product for model evaluation and introduce constraints in addition to traditional runoff data.

Noah-MP LSM
We employed the Noah-MP (Multi-Physics) LSM in this study.Noah-MP is an augment version of the Noah LSM [34] developed by Niu et al. [35] and Yang et al. [36] to improve the representations of various processes, including vegetation, snow [37], frozen soil [38] and hydrology schemes [39,40], etc.Additionally, Noah-MP has different parameterizations for several key land surface processes (Table 2); therefore, it can be used to conduct ensemble model simulations.Various studies have confirmed the performance of Noah-MP in reproducing energy and water fluxes at various scales and under different hydrometeorological conditions [41][42][43].
A snowpack can be divided into up to three layers in Noah-MP, and a series of processes including accumulation, compact, and melting are resolved by physically-based schemes (see [35] for details).More specifically, snowmelt is simulated by energy balance scheme and melt in the ith snow layer occurs when the internal temperature exceeds the freezing threshold temperature.The meltwater flows downward to the next snow layer and eventually reaches the land surface, where it is treated in the same way as rainfall by the model.
Runoff modeled by Noah-MP comprises both surface and subsurface runoff, which are represented mainly by a TOPMODEL-based saturation-excess runoff scheme [39].This scheme is suitable in the UCRB because of its humid and semi-humid climate.

Model Setup
The spatial and temporal resolutions of the modeling were determined based on the CMFD forcings, i.e., 0.1 • and 3 h.To match the spatial resolution of the CMFD data, we spatially resampled the land cover and soil datasets from their original resolutions to 0.1 • , using a majority resampling technique that assigns a dominant land cover/soil type to each 0.1 • grid cell.We downloaded the default parameters for each land cover/soil type from http://www.ral.ucar.edu/research/land/technology/noahmp_lsm.php.These parameters were calibrated globally using Modified IGBP MODIS 20-category vegetation (land-use) data and Hybrid STATSGO/FAO soil texture data.Options selected for key land surface processes parameterizations are listed in Table 2.
Before calibrating and validating the model, we spun up the model by repeatedly running 100 times simulations using 1979 and 1980 atmospheric forcings.This process created an equilibrium state of the model and reproduced a snow "reservoir" that defined as glaciers in the model [44], since glaciers are not treated in detail as it was done within the work of Prasch et al. [45], and Mauser and Prasch [46].

Model Calibration and Validation
In this paper, we mainly focus on the runoff and snow processes because they are directly related to our analysis of MR and MCR.
Note that surface runoff routing is not parametrized in the Noah-MP model; however, the bias is negligible because we mainly focus on runoff at the monthly scale and the mean residence time of surface runoff is less than 10 days.At each runoff gauging station, we identified the grid cells that comprise the associated drainage area and extracted the corresponding modeled runoff.The spatially averaged runoff in these grid cells is compared to the measured runoff for calibration.
Noah-MP is capable of modeling a variety of snow-related variables including the snow water equivalent (SWE), snow depth and SCF.However, few in-situ measurements are publicly available in the study area.Therefore, we adopted a remotely-sensed snow products for calibration.Comparing with SWE and snow depth products, the SCF products provides a relative accurate estimation [47,48].Therefore, we performed SCF calibration at the sub-region scale by spatially averaging the SCF values from Noah-MP and MOD10CM.
Due to the large computational demand and limited resources, it is impossible to perform automatic calibration in this study.Alternatively, based on two indices, i.e., the Nash-Sutcliffe coefficient (NSE) and correlation coefficient (r), we relied on the traditional trial and error approach for calibration.
where Q obs,i and Q sim,i are the observed and modeled runoff (mm•month −1 ) at time step i, respectively, and n is the total time step.Table 3 shows the calibrated parameters based on the observed runoff (from 1981 to 2000) and the MODIS derived SCF (from 2001 to 2005).According to the work of Cai et al. [41] and Niu et al. [39], melting factor (M) (see [37] for details) and threshold temperature of precipitation phases (T s ) are two sensitive parameters controlling snow dynamics while decay factor ( f ), maximum subsurface runoff coefficient (R sb,max ) and soil saturated hydraulic conductivity (K sat ) are sensitive parameters in runoff simulations.Due to the computational constraints, we only calibrated a lumped value (area average) in UCRB rather than spatially distributed values.We tested different multiplication factors of K sat , which has a default spatial distribution (see model setup for details), to determine the optimized values.
The performance of the model in the UCRB was further validated using runoff (from 2001 to 2010) and MODIS SCF (from 2006 to 2010) data.To provide a comprehensive evaluation, we also included the ET_Zhang product for comparison.Based on the parameterization schemes of different models, there are generally two approaches to quantifying MR and MCR.Models such as SRM [49] simulate the runoff generated by meltwater separately, from which MCR can be easily estimated [3,22].In contrast, models such as Noah-MP compute the melt and rainfall runoff together, making it difficult to separate MR from total runoff.However, models can be used to conduct numerical experiments and isolate behaviors that are difficult to isolate in nature.To quantify MCR in the UCRB, we conducted two numerical experiments using Noah-MP.We treated the validated simulation as a baseline experiment (Exp0) that the water flux incident at the soil surface consists of both rainfall and simulated meltwater.We then ran the model for a second time (Exp1) with simulated meltwater removed from the grid cell (we set meltwater to 0 after it was simulated by Noah-MP), and under such condition, the water flux incident at the soil surface is equal to rainfall.
We assume that the difference in runoff between the two experiments (Exp1-Exp0) is due to MR, and MCR can be estimated accordingly: where R Exp0,i and R Exp1,i are the runoff values (mm•month −1 ) simulated by Exp0 and Exp1, respectively, in the ith month.

Statistical and Trend Analyses
In this paper, we adopted two indices to evaluate the timing of MR and MCR.The first included the peak time (PT) for both MR and MCR, which is determined by identifying the month corresponding to the maximum MR or MCR in the year.The other index is "center time" (CT) for MR only.CT has been widely used to quantify the time of MR in related studies [17,50,51] as follows: where t i is time in months starting from January, q i is the corresponding melt runoff (mm/month) in the ith month, and n is the total number of months.
Although Stewart et al. [17] used total runoff from each water year to derive CT, we used MR directly because the UCRB is not dominated by snow and also because we can obtain MR values from numerical simulations, which can provide more accurate information regarding the timing.
With respect to trend detection, we employed the non-parametric Mann-Kendall (M-K) [52,53] test to evaluate the long-term trends in MR and MCR.Both Theil-Sen trend estimator and M-K trend significance (P) can be estimated for a given time series.In this paper, we only considered P > 0.95 as statistically significant.

Model Calibration and Validation
We first calibrated the snow processes by matching the MOD10CM SCF products in the calibration period from 2001 to 2005. Figure 2a-j (left part) shows the calibrated SCF values, and Table 3 shows the optimized parameters of M and T s .The figure shows that the model produces good results in terms of the correlation coefficient r with a maximum value of 0.85 for MJ, which indicates that the modeled SCF varies in accordance with MOD10CM SCF.However, the NSE coefficient reflected lower accuracy and the NSE values of JSJSY, JSJXY and CSH were negative.The relatively low NSE values can be attributed to both observation and model biases.Pu et al. [47] reported an overall accuracy of 90% for MODIS 8-day composite SCF estimation but a relatively high bias in regions where the snow depth was shallow, which is the case in JSJXY and CSH.In the validation period from 2006 to 2010 (Figure 2a-j (right part)), a good agreement between the modeled and observed values was further confirmed by both NSE and r, and no negative NSE avalues were produced during this period.To provide a clearer comparison, we plotted the climatology of SCF averaged over the entire 10 years beginning from 2001 (Figure 2k-t).These plots show that the modeled SCF can capture the overall variation in SCF estimated from MODIS data in all 10 sub-regions considering the challenge of modeling snow processes in the study area and the uncertainties associated with the remotely sensed SCF product.
After calibrating the snow-related parameters, we calibrated the runoff-related parameters using observed runoff data.Because the temporal coverage of the runoff data varied at six runoff gauging stations (Table 1), we chose the Gaochang, Beibei, and Yichang stations to perform the calibration using data from 1981 to 2000.Runoff data collected at all six stations beginning from 2001 were then used for validation.Figure 3c,e,f (left part) illustrate that the model can effectively reproduce the runoff variations during the calibration period at the Gaochang, Beibei and Yichang station in terms of both NSE (∼ 0.88) and r (∼ 0.95).The validation period (Figure 3a-f (right part)) exhibits good agreements between the modeled and observed values; thus, the parameters we calibrated are suitable in the study area (Table 3).Runoff can reflect the interactions between various hydrological processes, and accurate runoff simulations can indirectly reflect the performance of these uncalibrated processes.The NSE and r values are 0.91 and 0.98 at the Yichang station, which represents the integrated runoff from the entire UCRB.To provide a more comprehensive validation of the modeling results, we used the ET_Zhang dataset to validate the ET simulation.Figure 4a-i shows overall good agreements between modeled ET and the ET_Zhang estimated ET in ten sub-regions based on RMSE and r values.ET and runoff are two major components of precipitation partition, and the accurate simulations of ET and runoff reflect the effectiveness of the model we calibrated for the UCRB.

Overall Statistics of Melt Runoff
Figure 5a-j presents the monthly MCR calculated in the 10 sub-regions during the entire study period.The figure shows that JSJSY, YLJ, DDH and MJ are the four sub-regions with the most significant melt runoff.These four sub-regions are located in regions of the QTP with relatively high elevation and low temperatures.However, the MCR in CSH and WJ is negligible.In the following sections, we mainly focus on analyses of JSJSY, YLJ, DDH and MJ.
Based on Figure 5, we estimate the average MCR values in JSJSY, YLJ, DDH and MJ are 6.0%, 3.9%, 5.3% and 4.0%, respectively, from 1981 to 2010.On average, meltwater accounts for 2.0% of natural runoff at the Yichang runoff gauging station.

Temporal Characteristics of Melt Runoff
As shown in Figure 5a-j, the MCR values of different sub-regions shows significant inter-month variability, indicating the seasonal variation of MCR.We plotted the monthly climatology MCR values based on 30 years of simulation (Figure 5k-t) to illustrate the temporal pattern of MCR on an annual basis.In all sub-regions, MCR generally starts to increase in February (Figure 5k-p) or January (Figure 5q-t) which corresponds with the snow melt season in this area.As reported by Qin et al. [54], snow cover on the QTP normally begins to appear in mid-September and reaches an initial peak in January, which is in accordance with our SCF results shown in Figure 2k-n.In some sub-regions, two MCR peaks occur during a year in JSJSY, YLJ, DDH.The first peak occurs in July, June and May, respectively in these three sub-regions, with MCR values being as large as 15%.By further analyzing the corresponding climatology temporal variation in MR (Figure 5k-m; black dashed line), we found that the lag time between the MCR peaks of JSJSY and YLJ resulted from MR variation, which is related to the air temperature differences.However, the lag time between YLJ and DDH is mainly due to the variation in rainfall runoff since their MR peaks are both in June.In addition to the notable peak during the summer, a small winter peak can be observed in JSJSY, YLJ, and DDH around November.In the other seven sub-regions, a single MCR peak can be observed between March and May based on mean air temperature differences.
To investigate the long-term variation in MCR, we computed the annual MCR values in JSJSY, YLJ, DDH, and MJ (Figure 6a-d; blue solid line), and these values were used to perform M-K tests.The results show that all these four sub-regions experienced MCR increases during the study period (Figure 6a-d; red dashed line); however, only MJ shows a statistically significant increasing trend at 0.11%/year (P > 0.95).The PT of MCR shows no obvious variation at the monthly time scale (Figure 6e-h; red dashed line).
We also performed the M-K tests for MR.The long-term trend (Figure 7a-d; red dashed line) is different from that of MCR.Both JSJSY and MJ exhibits an statistically significant increasing trend of 0.02 mm/year and 0.04 mm/year, respectively.The PT of MR in all four sub-regions remains unchanged (Figure 7e-h; red dashed line), while CT of MR shows a slight and insignificant decreasing trend that reflects the advance of melt runoff (Figure 7i-l; red dashed line).

Spatial Characteristics of Melt Runoff
Figure 8 presents the monthly averaged spatial distribution of MCR, from which the melt runoff development can be easily identified.As previously shown in Figure 5, substantial MCR occurs between April and July.By comparing Figure 8 to Figure 1, we found that the MCR in the study area generally corresponded to the distribution of glaciers near Tanggula Mountain Range and Hengduan Mountain Range.
In JSJSY, MCR begins to increase in March in the downstream area near the Shigu runoff gauging station.It can be attribute to the relatively low elevation of this area, resulting in higher air temperatures which trigger the melting process earlier.The MCR in the upper reaches begins to increase in June, when most melting occurs in the glacierized area of the Tanggula Range (e.g., Mt.Geladandong).This MCR pattern is in accordance with the results reported by [23] who found that the continental-type glaciers in this region generally melt from June to August.Due to the high temperature in July, widespread melting occurs in JSJSY and causes the first MCR peak during the year.Then, October melting in the middle reaches near Mt.Shaluli leads to the second peak.
In YLJ, melting begins near the glacierized Mt.Gongga Region in March and in the upper reaches in May.MCR develops very quickly and peaks in June before decreasing.Then, MCR near two glacierized areas, Mt.Shaluli and the Mt.Gongga Region of Hengduan Mountain Range, starts to increase again in October, contributing to the second peak of MCR in this region.
In DDH and MJ, MCR generally occurs near the glacierized area of the Hengduan Mountain Range (e.g., Mt.Gongga and Mt.Yaomei) in April, July, and August.Additionally, widespread melting occurs in May and June, generally in the upper reaches.MCR values in these two sub-regions are very small after August.

Importance of Melt Runoff in the Upper Changjiang River Basin
We used a process-based LSM to quantify MCR in the UCRB for the period between 1981 to 2010.The overall multi-year average MCR in the entire UCRB is about 2%, mainly due to large downstream areas which receive little snow.However, sub-regions JSJSY, YLJ, DDH, and MJ have relatively large multi-year average MCR values of 6.0%, 3.9%, 5.3%, and 4.0%, respectively.Thirty years of monthly climatology MCR values suggest that the melt runoff constitutes up to 15% of the total runoff during summer, indicating the important role of melt runoff in these areas.The MCR becomes even more important during warm and dry years.For example, the MCR in JSJSY reached a notable peak value of 37.25% during 2001 when the precipitation was below the normal level [55].As such, meltwater plays an important role in regulating river runoff and sustaining a more stable runoff variation [56].Moreover, all major rivers in JSJSY, including the Togto, Dam Chu, Garchu, and Bi Chu, originate from and are supported by glaciers located in the source region [23].For this reason, MR helps sustain the well-being of the ecosystems and environment in this region.
Previous studies of melt runoff mainly focused on the source region (e.g., [20][21][22][23][24]), because it is regarded as the "water tower" [57].However, our results indicate that YLJ, DDH, and MJ have MCR values that are comparable to those of the source region.Moreover, based on our results that CT is advancing in the study area, MR should be taken into consideration for the flood control operations of these reservoirs to reduce the release of useful melt runoff during winter and early spring [17].

Effects of Climate and Land Cover Changes
As confirmed by many studies e.g., ([58-60]), climate change in the UCRB has already influenced the glaciers and associated runoff in this region [22,24,56].
Our study found a statistically significant increasing trend of both MR and MCR in MJ.Previous studies generally suggested that the increasing temperature and decreasing precipitation were responsible for the decreased runoff in MJ [61,62].A warmer climate tends to intensify the melting process of glacier/snow and increase MR.These factors, together with the precipitation decrease, increased MCR, which is consistent with our results.
Unlike the trends observed in MJ, our study indicates that MR increases and MCR decreases in JSJSY, which can be attributed to the increases in temperature and precipitation in this region [63], i.e., the increase in MR is offset by an increase in rainfall runoff.Additionally, glaciers in JSJSY are mainly cold-type or polar-type with relatively slower dynamic responses to global warming [24].
Additionally, we observed a slight advancement of CT, which is consistent with the overall warming in the UCRB.As reported by Wang et al. [51], shift can be observed in the snow melt period in the source region of the Changjiang River.The start time of melting advances 0.9-3 days/10a, or approximately 0.003-0.01months/year, which is in accordance with our results in JSJSY.However, the current trend is not statistically significant according to our results.
Besides climate, MCR is also influenced by land cover.On one hand, previous studies showed that vegetation affects snow accumulation by interception and sublimation, and affects snowmelt by alerting the energy balance of land surface [64][65][66].However, the relationship between snow and vegetation is highly variable and complicated, depending on climatic, topographic and vegetation conditions [65,67,68]; On the other hand, vegetation modifies runoff yield mainly by interception and ET [69,70].The effects of land cover change on MCR cannot be evaluated since we adopted the constant land cover during the entire study period.Moreover, from Section 3.4 we can see that glaciers and elevation have more significant effects on MCR.

Uncertainties in this Study
MR and MCR are difficult to measure at large scale, such as the UCRB scale and few studies have focused on the entire UCRB.Previous studies mainly focused on the source region of Changjiang.
For comparison, we also evaluated the MCR in the the source region and obtained a value of 8.9%, which is similar to those reported by Liu, et al. [22], Yang and Ming-Ko [71], and Shen et al. [23] of 11%, 9.2%, and 9%, respectively.However, other studies presented different results.For example, Zhang et al. [20] reported that MR accounted for up to 28.7% of the total runoff in the source region.Therefore, these results suggest that uncertainties may exist in quantifying MCR based on numerical models.

Atmospheric Forcings
Accurate atmospheric forcings are critical for land surface modeling [72], and precipitation is the most sensitive variable [73,74].The CMFD forcings we used in this paper have been widely used for modeling in China due to their high spatiotemporal resolution and relatively high accuracy.To evaluate the overall accuracy of the forcings in the UCRB, we compared the CMFD fields to corresponding measurements fom the China Meteorological Administration (CMA) ground observational sites.Eighty CMA sites are located in our study area, and data from these sites can be obtained online (http://data.cma.cn/data/detail/dataCode/SURF_CLI_CHN_MUL_DAY_CES_V3.0.html).Figure 9 illustrates the grid-to-site comparison of air temperature and precipitation from 1981 to 2010 based on extracting the CMFD values according to the locations of CMA sites at both the daily and monthly time scales.The results show that although the CMFD air temperatures agree well with those of CMA sites, the CMFD precipitation is biased, especially at the daily time scale.However, by comparing with the study evaluating the Global Land Data Assimilation System (GLDAS) forcings over China [75], we found CMFD forcings providing better air temperature and precipitation fields over the study area, in both accuracy and resolution of the data.Additionally, our analyses in this paper are performed at monthly or annual time scales, therefore the overall accuracy of analyses can still be guaranteed.

Snow-Related Parameterizations
Although the overall parameterizations of snow in Noah-MP are physically based, several simplifications and generalizations of processes still exist in the model [35,37], which is the case in many other LSMs [76].Studies have documented the wide spread of model results [19,76,77] when applied to the same numerical experiments, indicating the discrepancies and uncertainties in different parameterizations implemented.Regarding MCR, Siderius et al. [15] reported different MCR of Ganges estimated from four different models, particularly in upstream area.Essery et al. [77] suggested that proper model calibration can benefit the snow simulations in many models.In this study, we calibrated two most sensitive snow-related parameters using the remotely-sensed MODIS SCF data.However, due to the lack of ground-based snow data (e.g., snow water equivalent and snow depth), it was impossible to calibrate other parameters e.g., liquid water holding capacity, which are important to snow processes [37].

Quantification of Melt Water Contribution to Total Runoff
In this paper, we conducted two numerical experiments to isolate the MR and estimate MCR.The actual hydrological processes are highly nonlinear and complex [78]; therefore, the direct removal of meltwater from the simulations may affect processes such as infiltration and ET and contribute to the runoff differences between the two experiments.However, we assume the bias in this approach is small in the UCRB because the meltwater volume is relatively small compared to that of rainfall runoff.
Studies have proposed different approaches to quantify MCR; however, uncertainties also exist in these approaches.For example, Immerzeel et al. [3] adopted the Normalized Melt Index (NMI) to investigate MCR, but the runoff coefficient of MR was difficult to estimate because observed data were generally unavailable.Additionally, Siderius et al. [15] used the ratio of melt inputs to rainfall inputs as an indicator, but meltwater of previous month may have been incorrectly considered as rainfall in the subsequent month.

Conclusions
In this paper, we adopted the Noah-MP LSM to investigate the spatiotemporal characteristics of MR and MCR in UCRB.Remotely-sensed SCF (2001-2005) and ground-measured runoff (1981-2001, whenever available) data were used to calibrate five most sensitive parameters in the Noah-MP model.Further validation against SCF, runoff and ET observations confirmed the performance of the Noah-MP model in the study area.We then used the calibrated model to conduct two numerical experiments to quantify MR and MCR in UCRB from 1981 to 2010, and the following conclusions were drawn from this study.1.On average, MR accounted for about 2.0% of natural runoff at the Yichang gauging station, which controls the entire UCRB.More specifically, JSJSY, YLJ, DDH and MJ, the four sub-regions of UCRB, had the most significant MCR values of 6.0%, 3.9%, 5.3% and 4.0%, respectively during the study period.By comparison, the MCR of values in CSH and WJ were negligible.2. The MCR values in different sub-regions showed significant seasonal variability.From a monthly climatology perspective, MCR generally begins to increase in January or February in all sub-regions.In JSJSY, YLJ, and DDH, double MCR peaks were observed in summer and winter while only a single peak was observed in the remaining sub-regions.The PT varied by sub-regions, mainly due to air temperature gradients.3. M-K tests indicated that the annual MCR values of JSJSY, YLJ, DDH, and MJ showed increasing trends during the study period; however, only that of MJ was statistically significant, with a value of 0.11%/year (P = 1.00).The annual MR values of both JSJSY and MJ displayed statistically significant increasing trends of 0.02 mm/year (P = 0.96) and 0.04 mm/year (P = 0.98), respectively.No trends were found for PT values of MR and MCR, in contrast, advancing trend were observed for th CT values of MR in JSJSY, YLJ, DDH, and MJ, and the corresponding values were 0.01 months/year, 0.02 months/year, 0.01 months/year and 0.01 months/year, respectively.4. The spatial distribution of MCR values in the study area generally consistent with the distribution of glaciers and elevation.In JSJSY, YLJ, DDH, and MJ, MCR generally begins to increase near glacierized areas of Tanggula Range and Hengduan Mountain Range.MCR was then influenced by the melting of the seasonal snow pack in the upper stream area. 5.The uncertainties in this study may result from the atmospheric forcings, the snow related parameterizations in the Noah-MP model, and the approach we adopted to quantify the MCR.
According to our study, MR is important in the UCRB, especially in JSJSY, YLJ, DDH, and MJ, which are located on the QTP, especially during warm and dry years.Climate change has already affected the snow-related processes in the UCRB, especially in MJ, which exhibited significant increase in both MR and MCR.Melt water should be considered in flood control operations and water resource management in these regions.

Figure 1 .
Figure 1.Map of the Upper Changjiang River Basin (see Section 2.1 for details).Blue shaded areas are the glaciers extracted from "The second Chinese glacier inventory" [27].

Figure 2 .
Figure 2. Comparison between the modeled and remote sensing retrieved (MOD10CM) spatial average snow cover fraction (SCF) in then sub-regions (a-j) during the calibration period (2001 to 2005) and the validation period (2006 to 2010).The monthly climatology SCF averaged between 2001 and 2010 (k-t) is also shown.

Figure 3 .
Figure 3.Comparison between modeled and observed runoff at six runoff gauging stations (a-f) during the calibration period (1981 to 2000) and validation period (2001 to 2010).

Figure 5 .
Figure 5. Modeled monthly melt runoff contribution to total runoff (MCR, %) in ten sub-regions (a-j, blue line) from 1981 and 2010.The corresponding monthly climatology MCR values averaged over the same period (k-t, blue line) are also shown, and the corresponding monthly climatology melt runoff (MR, mm/month) is illustrated for comparison (k-t, black line).

Figure 6 .
Figure 6.Modeled interannual variation in (a-d, blue solid line) melt runoff contribution to total runoff (MCR, %) and (e-h, blue solid line) the corresponding peak time of MCR (PT, month) in four sub-regions (JSJSY, YLJ, DDH, and MJ) with the most notable MCR values from 1981 to 2010.The Mann-Kendall (M-K) trends for MCR (a-d, red dashed line) and PT (e-h, red dashed line) are also shown, and the Theil-Sen trend estimator (Trend) and corresponding M-K trend significance (P) are listed above each panel.

Figure 7 .
Figure 7. Modeled interannual variation in (a-d, blue solid line) melt runoff (MR, mm/year), (e-h, blue solid line) the peak time of MR (PT, month), and (i-l, blue solid line) the center of MCR (CT, month) in the 4 sub-regions (JSJSY, YLJ, DDH, and MJ) with the most notable MCR values from 1981 to 2010.The Mann-Kendall (M-K) trend calculated for MR (a-d, red dashed line), PT (e-h, red dashed line), and CT (i-l, red dashed line) are also shown, and the Theil-Sen trend estimator (Trend) and coresponding M-K trend significance (P) are listed above each panel.

Figure 8 .
Figure 8. Modeled monthly spatial distribution of melt runoff contribution to total runoff (MCR, %) in the entire Upper Changjiang River Basin from 1981 to 2010.

Figure 9 .
Figure 9. Scatter plot between China Meteorological Administration (CMA) ground sites measured and the corresponding China Meteorological Forcing Dataset (CMFD) grid retrieved variables of (a) daily and (b) monthly air temperature; and (c) daily and (d) monthly precipitation.

Table 1 .
[27]acteristics of the ten sub-regions in the Upper Changjiang River Basin.Dominant land cover type was retrieved from[28], glacier area was retrieved from[27], and temporal series of runoff are shown at six runoff gauging stations.

Table 2 .
Options selected for key land surface processes parameterizations in Noah-MP.

Table 3 .
The parameters calibrated in this study for the snow cover fraction and runoff, and the optimized parameters values.