Modeling Ecohydrological Processes and Spatial Patterns in the Upper Heihe Basin in China

The Heihe River is the second largest inland basin in China; runoff in the upper reach greatly affects the socio-economic development in the downstream area. The relationship between spatial vegetation patterns and catchment hydrological processes in the upper Heihe basin has remained unclear to date. In this study, a distributed ecohydrological model is developed to simulate the hydrological processes with vegetation dynamics in the upper Heihe basin. The model is validated by hydrological observations at three locations and soil moisture observations at a watershed scale. Based on the simulated results, the basin water balance characteristics and their relationship with the vegetation patterns are analyzed. The mean annual precipitation and runoff increase with the elevation in a similar pattern. Spatial patterns of the actual evapotranspiration is mainly controlled by the precipitation and air temperature. At the same time, vegetation distribution enhances the spatial variability of the actual evapotranspiration. The highest actual evapotranspiration is around elevations of 3000–3600 m, where shrub and alpine meadow are the two dominant vegetation types. The results show the mutual interaction between vegetation dynamics and hydrological processes. Alpine sparse vegetation and alpine meadow dominate the high-altitude regions, which contribute most to the river runoff, and forests and shrub contribute relatively small amounts of water yield.


Introduction
Vegetation, topography and climate variables interactively influence multi-scale ecohydrological processes in a large catchment.Vegetation patterns are a key factor affecting the water balance and catchment water yield.Previous studies have demonstrated the impact of forest change on river discharge in China [1,2].However, few studies have addressed the relationship between hydrological processes and vegetation patterns.In particular, the relationship between forest and water in the headwater catchments of arid inland basins remains unclear [3].In arid inland basins, the river discharge generated from the mountain regions greatly affects the socio-economic development and ecosystem sustainability of the downstream regions.Therefore, understanding the complex relationship between spatial vegetation patterns and hydrological processes is important for integrated river basin management.Hydrological observations and small-scale experiments provide limited information for understanding the spatial patterns of ecohydrological processes in a large river basin [4].
Instead, a physically based distributed model that links ecohydrological processes across scales is needed to analyze the spatial variability of hydrology with vegetation pattern.
In the past few decades, with global climate change and population growth, the water shortage and ecosystem degradation in many river basins have gained increasing concern worldwide.Meanwhile, changes in natural river runoff in the headwater catchments worldwide have been reported in previous studies [5][6][7][8][9][10][11].Because of the less direct influence of human activities in headwater areas, runoff changes are mainly caused by the mutual interactions among climate, vegetation and hydrology.Yang et al. [12] reported that climate contributions to river runoff showed a large spatial variation over China, and several studies attempted to attribute runoff change to the impacts of climate and vegetation changes [13][14][15].However, vegetation is absent or simply parameterized in traditional hydrological models [16].It is relatively difficult to evaluate the influence of vegetation change on runoff due to data availability and methodological limitations [17,18].Therefore, it is important to develop the ecohydrological models for understanding and predicting changes in the regional water availability under the changing environment.
Modeling hydrological processes at the catchment scale requires a flexible distributed scheme to represent the catchment topography, river network and vegetation patterns [19].Previous studies of distributed hydrological models focused on the representation of the heterogeneity of catchment landscape and physically based descriptions of hydrological processes, such as the Systeme Hydrologique Europeen (SHE) and the Variable Infiltration Capacity (VIC) model [20,21].Yang et al. [22,23] proposed the geomorphology-based hydrological model (GBHM) considering sub-grid parameterization, which has been employed in macroscale studies [24,25] and mesoscale studies [13,14,26,27].To better simulate the hydrological changes in a catchment, it is necessary to couple the vegetation pattern and vegetation dynamics in hydrological models [19].Land surface models, such as the second version of the Simple Biosphere model (SiB2) [28], have mainly focused on the role of vegetation in the water-heat-carbon cycle.However, the catchment hydrology was poorly represented in SiB2.Further research embedded the SiB2 into the GBHM and developed a distributed biosphere hydrological model called the Water and Energy Budget-Based Distributed Hydrological Model (WEB-DHM) [29].A comparative analysis of the WEB-DHM and GBHM has been carried out and indicated that the WEB-DHM could be used to predict streamflow and conduct the water and energy flux estimations [30].Several land surface models can simulate the cryosphere processes, such as the Community Land Model (CLM) [31].However, they typically simplify the catchment hydrological processes and cannot accurately simulate streamflow.Several studies attempt to develop physically based models to couple hydrological processes with ecological processes.Tague and Band [32] developed the Regional Hydro-Ecologic Simulation System (RHESSys) to simulate both hydrological processes and vegetation dynamics.However, in RHESSys, the frozen soil process is not considered and the soil temperature is estimated based on the average daily temperature.Maneta and Silverman [33] and Ivanov et al. [34] developed spatial distributed models to simulate the water, energy balance and vegetation dynamics.However, the frozen soil is not adequately considered, as there are only two soil thermal layers in this model.Therefore, it is important to develop a distributed model that couples hydrological processes, cryosphere processes and vegetation dynamics.
There are several previous studies focusing on the hydrological modeling in the Heihe River basin.Jia et al. [35] applied the water and energy transfer process (WEP) model to simulate and predict the annual runoff changes due to climate and land use changes.Wang et al. [36] and Zhang et al. [37] applied heat and water coupling models to a small experimental catchment in the upper Heihe basin to test the model applicability in the high, cold mountainous area.Zhou et al. [38] chose several modules for the hydrological processes in cold regions and linked them to a catchment hydrological model to improve the hydrological modeling capability in the cold mountainous area.Zang and Liu [39] applied the Soil and Water Assessment Tool (SWAT) in the Heihe basin to analyze the flow trends of green and blue water.Yang et al. [19] developed a distributed scheme for ecohydrological modeling in the upper Heihe River.However, this scheme needs further coupling of vegetation dynamics into the hydrological simulation and improvement of the cryosphere hydrological processes.Because of the complexity of ecohydrological processes in the upper Heihe basin, further research to develop a sophisticated ecohydrological model is still ongoing.A major research plan entitled "Integrated research on the ecohydrological processes of the Heihe basin" has been launched by the National Natural Science Foundation of China since 2010 [40].One of the integrated research projects in this research plan aims to develop a distributed ecohydrological model for the cold mountainous regions of inland basins.With the help of the Heihe Watershed Allied Telemetry Experimental Research (HiWATER), which is a comprehensive ecohydrological experiment in this research plan [4], there are new opportunities for ecohydrological model development in the Heihe basin.
The major objectives of this study are to (1) develop a distributed ecohydrological model by coupling hydrological processes with vegetation dynamics and ( 2) analyze the relationship between the water balance and vegetation patterns in the upper Heihe basin.
In the following sections, the features of the study area are presented, followed by the data descriptions.The distributed ecohydrological model is then introduced.In the results and discussion section, the model validation is presented, followed by the analysis of water balance characteristics and the spatiotemporal variability of runoff.Finally, the catchment ecohydrological pattern is analyzed from the water balance and its relation to the vegetation distribution.

The Upper Heihe Basin
Heihe basin is the second largest inland basin in Northwest China.The upper reach of the Heihe River, which is gauged at the Yingluoxia hydrological station (see Figure 1), has a drainage area of 10,005 km 2 .The upper Heihe basin generates nearly 70% of the total river runoff, which supplies agricultural irrigation and benefits the social economy development in the middle and lower basin [19,41].There are two tributaries in the upper Heihe basin, namely, the West Tributary and East Tributary, both of which originate in the Qilian Mountains.The annual precipitation is between 200 mm and 700 mm with highly seasonal variability, and nearly 60% of the total annual rainfall is concentrated in summer from June to September.The upper Heihe basin has an elevation of 1700-5200 m, low temperature and relatively abundant vegetation types [41][42][43][44].The major vegetation types include coniferous forest (Picea crassifolia Kom.), shrub (Potentilla fruticosa Linn.), steppe (Stipa purpurea Griseb), alpine meadow (Kobresia pygmaea Clarke), alpine sparse vegetation (Saussurea medusa Maxim.), and desert (Sympegma regelii Bunge) (see Figure 2).

Data Used in the Study
Two types of data are used in this study: the first type is the data used to build and run the ecohydrological model, which include the geographic information and climatic forcing data, and the second type is the data used for model validation.
Meteorological observations are available at several stations within the study basin and its surrounding areas, as shown in Figure 1.The data are acquired from the National Meteorological Information Center affiliated with the China Meteorological Administration [45].The observed meteorological data include daily precipitation, temperature, sunshine hour, wind speed and relative humidity.In this study, daily gridded precipitation is interpolated from the gauge data using the method developed by Shen and Xiong [46].Other forcing data are interpolated using the inverse distance method.The hourly temperature is estimated from the daily maximum and minimum temperature using a sine curve.Hourly precipitation is estimated from the daily data according to the duration.Duration of precipitation within a day is determined by the precipitation amount estimated from the regional historical records.The starting time of the precipitation in a day is decided randomly and the hourly precipitation is specified using a normal distribution.The wind speed and relative humidity are assumed uniform in a day.The location of glaciers and their areas and volumes are obtained from the first and second glacier inventory datasets of China [47][48][49], which is downloaded from the Cold and Arid Regions Science Data Center at Lanzhou [50].Digital elevation data with a resolution of 90 m resolution are downloaded from the Shuttle Radar Topography Mission (SRTM) database [51].The vegetation map (see Figure 2) is obtained from the Institute of Botany, Chinese Academy of Sciences [52,53].The leaf area index (LAI) is estimated from the remote sensing data by Fan [54].The soil map of Heihe basin is produced by the second national soil survey [55].We use the van Genuchten model to represent the retention curve of soil water.The soil water parameters used in this study, including the saturated hydraulic conductivity, residual soil moisture, saturated soil moisture, and parameters α and n in the van Genuchten model, are provided by the China soil dataset developed by Dai et al. [56].The spatial resolution of this dataset is 1 km.elevation data with a resolution of 90 m resolution are downloaded from the Shuttle Radar Topography Mission (SRTM) database [51].The vegetation map (see Figure 2) is obtained from the Institute of Botany, Chinese Academy of Sciences [52,53].The leaf area index (LAI) is estimated from the remote sensing data by Fan [54].The soil map of Heihe basin is produced by the second national soil survey [55].We use the van Genuchten model to represent the retention curve of soil water.The soil water parameters used in this study, including the saturated hydraulic conductivity, residual soil moisture, saturated soil moisture, and parameters α and n in the van Genuchten model, are provided by the China soil dataset developed by Dai et al. [56].The spatial resolution of this dataset is 1 km.To validate the runoff simulation, the stream discharge data are obtained from the Hydrology and Water Resources Bureau of Gansu Province.Soil moisture measured by the wireless sensor network (WSN) in the East Tributary is obtained from the HiWATER [57,58].The actual evapotranspiration dataset estimated by remote sensing is available in the study basin at a 1 km spatial resolution and monthly temporal scale [59,60], which is used for comparison with the simulated result.To validate the runoff simulation, the stream discharge data are obtained from the Hydrology and Water Resources Bureau of Gansu Province.Soil moisture measured by the wireless sensor network (WSN) in the East Tributary is obtained from the HiWATER [57,58].The actual evapotranspiration dataset estimated by remote sensing is available in the study basin at a 1 km spatial resolution and monthly temporal scale [59,60], which is used for comparison with the simulated result.

Representation of the Landscape
The distributed scheme used to develop the ecohydrological model is originally from the geomorphology-based hydrological model (GBHM) [61-64] and was improved by Yang et al. [19] for ecohydrological modeling.The distributed ecohydrological model developed in this study is called the geomorphology-based ecohydrological model (hereafter referred to as the GBEHM).The major development in the GBEHM includes replacing the evapotranspiration estimation from the Penman-Monteith equation to a simple biosphere model used in SiB2 [65] and adding an energy balance based module for simulating the cryosphere hydrological processes.
A grid system with a resolution of 1 km × 1 km is used to discretize the study catchment, and the river network is extracted from the 1-km digital elevation model (DEM), by which the study catchment is divided into 461 sub-catchments.The main rivers of 461 sub-catchments are used to represent the streamflow pathway of the study catchment, and the Horton-Strahler ordering system is used to define the flow routing sequence.The mean terrain properties (slope length and gradient) and a particular soil type are assigned to each 1-km grid and the area ratios of vegetation types are determined for each grid.
Furthermore, each 1-km grid cell is represented by a number of topographically similar "hillslope-valley" systems.Length, gradient and aspect of the hillslope are estimated from the 90-m DEM and averaged on the 1-km grid [19].The hillslopes within a 1-km grid are grouped according to vegetation type.The vertical structure of a hillslope is then subdivided into vegetation canopy, soil and bedrock.

Simulation of Ecohydrological Processes
Hillslope is the basic unit of ecohydrological simulation in the GBEHM, in which vegetation dynamics are coupled with the hillslope hydrological processes.The descriptions of hydrological processes, which mainly include the transfer of energy, water and carbon dioxide in the soil-plantatmosphere continuum, are especially designed for better simulation in the cryosphere.
During the growing season, the photosynthesis process of vegetation is simulated together with canopy energy transfer and canopy evapotranspiration.The radiation transfer in the vegetation canopy layer includes interception, reflection, transmission and absorption and is described by the same scheme in SiB2 [28,66].The energy balance equation in the canopy layer is expressed as [28]

Representation of the Landscape
The distributed scheme used to develop the ecohydrological model is originally from the geomorphology-based hydrological model (GBHM) [61-64] and was improved by Yang et al. [19] for ecohydrological modeling.The distributed ecohydrological model developed in this study is called the geomorphology-based ecohydrological model (hereafter referred to as the GBEHM).The major development in the GBEHM includes replacing the evapotranspiration estimation from the Penman-Monteith equation to a simple biosphere model used in SiB2 [65] and adding an energy balance based module for simulating the cryosphere hydrological processes.
A grid system with a resolution of 1 km ˆ1 km is used to discretize the study catchment, and the river network is extracted from the 1-km digital elevation model (DEM), by which the study catchment is divided into 461 sub-catchments.The main rivers of 461 sub-catchments are used to represent the streamflow pathway of the study catchment, and the Horton-Strahler ordering system is used to define the flow routing sequence.The mean terrain properties (slope length and gradient) and a particular soil type are assigned to each 1-km grid and the area ratios of vegetation types are determined for each grid.
Furthermore, each 1-km grid cell is represented by a number of topographically similar "hillslope-valley" systems.Length, gradient and aspect of the hillslope are estimated from the 90-m DEM and averaged on the 1-km grid [19].The hillslopes within a 1-km grid are grouped according to vegetation type.The vertical structure of a hillslope is then subdivided into vegetation canopy, soil and bedrock.

Simulation of Ecohydrological Processes
Hillslope is the basic unit of ecohydrological simulation in the GBEHM, in which vegetation dynamics are coupled with the hillslope hydrological processes.The descriptions of hydrological processes, which mainly include the transfer of energy, water and carbon dioxide in the soil-plant-atmosphere continuum, are especially designed for better simulation in the cryosphere.
During the growing season, the photosynthesis process of vegetation is simulated together with canopy energy transfer and canopy evapotranspiration.The radiation transfer in the vegetation canopy layer includes interception, reflection, transmission and absorption and is described by the same scheme in SiB2 [28,66].The energy balance equation in the canopy layer is expressed as [28] where C is the effective heat capacity of the canopy (J¨m ´2¨K ´1), T c is the canopy temperature (K), t is the time (s), R n is the absorbed net radiation (W¨m ´2), H is the sensible heat flux (W¨m ´2), λ is the latent heat of vaporization (J¨kg ´1), E is the evapotranspiration rate (kg¨m ´2¨s ´1), and ξ is the energy transfer due to water phase change (W¨m ´2), which is caused by the melting/freezing of the snow intercepted by the canopy.Thus, the energy budget is linked with water balance by the evapotranspiration component.The absorbed net radiation is calculated using the radiation transfer model described by Sellers et al. [65].The sensible heat flux of the canopy H is calculated as where ρ is the air density (kg¨m ´3), c p is the specific heat of air (J¨K ´1¨kg ´1), T c is the canopy temperature (K), T a is the air temperature (K), and r b is the bulk canopy boundary layer resistance (s¨m ´1).
The actual evapotranspiration rate of the vegetation layer is influenced by the leaf stomatal conductance, which is related to environmental properties, such as the canopy temperature, carbon dioxide concentration and water vapor in the air.In this study, canopy conductance is estimated as [28] where g c is the canopy conductance (m¨s ´1), m and b are the empirical coefficients of vegetation types, C s is the carbon dioxide partial pressure at the leaf surface (Pa), h s is the relative humidity, p is the atmospheric pressure (Pa), L T is the total leaf area index (m 2 ¨m´2 ), and A c is the canopy photosynthesis rate (mol¨m ´2¨s ´1) and is calculated as [67] A c " A n0 Π (4) where A n0 is the net assimilation rate for the leaves at the top of the canopy (mol¨m ´2¨s ´1), FPAR is the fractional interception of photosynthetically active radiation, and k is the time-mean extinction coefficient for photosynthetically active radiation.The net assimilation rate of the leaf is calculated as where A n is the leaf net assimilation rate (mol¨m ´2¨s ´1), A is the leaf photosynthetic rate (mol¨m ´2¨s ´1), and R d is the leaf respiration rate (mol¨m ´2¨s ´1).A and R d are estimated by the functions of the maximum catalytic capacity of the photosynthetic enzyme (V max ), canopy temperature and other environmental factors [28,68,69].Details about the calculation of g c , A and R d were given by Sellers et al. [28].Based on the canopy conductance g c , the canopy transpiration rate is expressed as [28] λE c " where E c is the canopy transpiration rate (kg¨m ´2¨s ´1), e* is the saturation vapor pressure (Pa), e a is the vapor pressure in the canopy (Pa), r b is the canopy boundary layer resistance (s¨m ´1), ρ is the density of air (kg¨m ´3), c p is the specific heat of air (J¨kg ´1¨K ´1), and γ is the psychrometric constant (Pa¨K ´1), W c is the canopy wetness-snow cover fraction, and the other parameters are the same as in Equations ( 1) and ( 3).The soil evaporation rate E g of the surface soil layer and the evaporation rate of the canopy interception E i are calculated using the same methods in SiB2 [28].E g and E i are added to E c as the total actual evapotranspiration over the hillslope unit.
In the non-growing season, ecohydrological simulation focuses on the cryosphere hydrological processes, which mainly include the snow melting process, soil freezing and thawing process, and glacier melting process.We use a similar scheme as in the Community Land Model version 4.0 (i.e., CLM4.0) [31] to represent the heat transfer in snow and frozen soil: where c is the volumetric snow/soil heat capacity (J¨m ´3¨K ´1), T s is the temperature (K) of the snow/soil layers, z is the vertical depth of snow/soil (m), and K T is the thermal conductivity (W¨m ´1¨K ´1).Equation ( 8) solves the snow/soil temperature with the boundary condition as the heat flux into the top surface layer of the snow/soil and zero heat flux at the bottom of the soil column.
The surface layer heat flux from the atmosphere is expressed as where h is the upper boundary heat flux into the snow/soil layer (W¨m ´2), S g is the solar radiation absorbed by the top layer (W¨m ´2), L g is the long-wave radiation absorbed by the ground (W¨m ´2), H g is the sensible heat flux from the ground (W¨m ´2), and λE g is the latent heat flux from the ground (W¨m ´2).
After solving Equation ( 8), the soil or snow layer temperature is evaluated to determine whether the phase change will take place [31].For the soil or snow layers, melting takes place under the condition of For the snow layers, freezing takes place under the condition of T i ă T f and W liq ą 0 (11) and for the soil layers, freezing takes place under the condition of In Equations ( 10)-( 12), T i is the temperature of the soil or snow layers (K), T f is the freezing temperature of water (K), W ice is the mass of ice in the soil or snow layers (kg¨m ´2), W liq is the mass of liquid water in soil or snow layers (kg¨m ´2), and W liq,max is the maximum mass of supercooled soil water, which is the liquid water that coexists with the ice (kg¨m ´2).The phase change rate is determined by the energy excess (or deficit) needed to change the temperature of soil or snow layers (T i ) to the freezing temperature (T f ).If the melting criteria in Equation ( 10) is met and energy excess is greater than zero, then the ice mass is calculated by [31] where W n`1 ice is the ice mass after the phase change (kg¨m ´2), W n ice is the ice mass before the phase change (kg¨m ´2), n + 1 and n refer to the time steps, U i is the energy excess (W¨m ´2), ∆t is the length of time step (s), and L f is the latent heat of ice fusion (J¨kg ´1).If the freezing criteria in Equations (11) or ( 12) is met and energy excess is less than zero, then the ice mass in the snow layers is adjusted by The ice mass in the soil layers is adjusted by and the mass of liquid water is adjusted by where W n`1 liq and W n liq are the mass of liquid water after and before the phase change (kg¨m ´2), respectively, in Equations ( 14)-( 16), and the other parameters are the same as in Equation ( 13).Moreover, glacier melting is simulated using the degree-day model [70,71].
The surface runoff is from the infiltration excess and saturation excess calculated by solving Richards' equation using an implicit finite difference method.The surface runoff flows through the hillslope into the stream via kinematic wave.The groundwater aquifer is treated as an individual storage corresponding to each grid.The exchange between the groundwater and river water is considered as steady flow and calculated using Darcy's law [64].The runoff generated from the grid is the lateral inflow into the river of the sub-catchment.The flow routing in the river network of the whole study catchment is solved using the kinematic wave approach: where Q is the discharge (m 3 ¨s´1 ), t is the time (s), x is distance along the river (m), A is the area of the cross-section (m 2 ), q is the lateral inflow to the river from the hillslope (m 3 ¨s´1 ), S 0 is the slope of the river bed, n r is the roughness of the river bed, and p is the wetting perimeter of cross-section (m).Equation ( 17) is solved using a nonlinear explicit finite difference method and Newton's iteration scheme.The time step of the GBEHM model is one hour.

Model Calibration and Performance Evaluation Metrics
Most model parameters are estimated from field observations or remote sensing data, and some vegetation parameters are specified from previous studies [19].The major parameters related to vegetation type are listed in Table 1.Considering the equilibration time for the hydrological state variables (e.g., soil moisture, soil temperature, groundwater table), a warm-up run of the model in the period of 1999-2000 is used to update the model state variables.We perform our analysis in the period from 2001 to 2012.The Nash-Sutcliffe efficiency (NSE) coefficient and relative error (RE) are used to evaluate the model performance.

Model Validation
The model has been validated using the soil moisture measured in the East Tributary of Heihe River, streamflow discharge observed at three hydrological stations (Figure 1) and actual evapotranspiration estimated from the remote sensing observations.At the catchment scale, the model is validated using the streamflow and spatial distribution of observed soil moisture.Figure 3 illustrates the monthly river discharge of the Yingluoxia, Zhamashike and Qilian hydrological stations, which are located at the outlet of upper Heihe River, the West Tributary and the East Tributary, respectively.The simulated streamflow generally shows a good agreement with the observed one.The NSE values for the three observation stations are 0.77, 0.80 and 0.67, respectively, and the absolute values of RE are smaller than 10%.The distributed soil moisture simulation of the top 5-cm layer is compared with the WSN observation at 4 cm from the surface of the HiWATER in the East Tributary of the upper Heihe basin.A comparison of the simulations and observations for the monthly averaged soil moisture over the East Tributary is shown in Figure 4.The areal average values of the simulated and observed soil moistures are highly similar at 0.41 and 0.40 in August, 0.40 and 0.39 in September, and 0.31 and 0.29 in October, respectively.The simulated soil moisture generally captures the spatial pattern of the observed soil moisture.However, the details need to be checked carefully in future research.
In addition, the simulated actual evapotranspiration (ET) is compared with the remote sensing-based estimation.Figure 5 illustrates the comparison of annual average ET between the model simulation and remote sensing-based estimation in the study catchment in the 2001-2012 period [59,60].Both of the estimations have close long-term basin average ETs and similar spatial patterns over the study catchment.The areal average ET is 310.8 mm for the model simulation and 306.7 mm for the remote sensing-based estimation.However, the remote sensing-based estimation shows a higher spatial variability than the GBEHM simulation.The remote sensing-based estimation mainly considers the energy balance, and the land surface roughness plays an important role.In addition, the ecohydrological model considers both energy and water balances when estimating the actual evapotranspiration, in which vegetation plays an important role.Vegetation may re-distribute the precipitation and moderate the spatial variability of soil moisture and evapotranspiration.

Water Balance Characteristics and Spatio-Temporal Variability of Runoff
Based on the simulation, the annual average water balance during 2001-2012 is calculated for the East and West Tributaries and the entire catchment, and the results are given in Table 2.For the upper Heihe basin, the annual average precipitation, actual evapotranspiration (ET) and runoff are 479.9mm, 310.8 mm and 169.0 mm, respectively.Comparing the water balance in the East and West Tributaries, the East Tributary has higher actual ET and runoff due to its higher precipitation.However, the West Tributary has a larger runoff ratio compared to the East Tributary due to the higher altitudes and relatively larger areas of glaciers.The seasonal characteristics of water balance for the East and West Tributaries and the entire catchment are also estimated and are shown in Table 3.In general, the East and West Tributaries and the entire catchment show similar seasonal patterns.For the entire catchment, the precipitation in winter (from December to February) is only 11.5 mm, and the actual evapotranspiration (ET) and runoff is also rather small (6.2 mm and 8.9 mm, respectively).The precipitation in spring (from March to May) is 81.4 mm, which is lower than the sum of the actual ET (66.0 mm) and runoff (19.3 mm).This implies that the runoff in spring is generated from not only precipitation but also snow and glacier melting.This characteristic is more obvious for the West Tributary due to the relatively larger area of glaciers (see Figure 2).As the glacier area is relatively small (less than 1% of the basin area), the contribution of glacier melting to river discharge in summer and autumn is quite low.In summer (from June to August), the precipitation is 281.6 mm, which is larger than the total of the actual ET (180.1 mm) and runoff (78.3 mm).This result implies that precipitation recharges the soil water and

Water Balance Characteristics and Spatio-Temporal Variability of Runoff
Based on the simulation, the annual average water balance during 2001-2012 is calculated for the East and West Tributaries and the entire catchment, and the results are given in Table 2.For the upper Heihe basin, the annual average precipitation, actual evapotranspiration (ET) and runoff are 479.9mm, 310.8 mm and 169.0 mm, respectively.Comparing the water balance in the East and West Tributaries, the East Tributary has higher actual ET and runoff due to its higher precipitation.However, the West Tributary has a larger runoff ratio compared to the East Tributary due to the higher altitudes and relatively larger areas of glaciers.The seasonal characteristics of water balance for the East and West Tributaries and the entire catchment are also estimated and are shown in Table 3.In general, the East and West Tributaries and the entire catchment show similar seasonal patterns.For the entire catchment, the precipitation in winter (from December to February) is only 11.5 mm, and the actual evapotranspiration (ET) and runoff is also rather small (6.2 mm and 8.9 mm, respectively).The precipitation in spring (from March to May) is 81.4 mm, which is lower than the sum of the actual ET (66.0 mm) and runoff (19.3 mm).This implies that the runoff in spring is generated from not only precipitation but also snow and glacier melting.This characteristic is more obvious for the West Tributary due to the relatively larger area of glaciers (see Figure 2).As the glacier area is relatively small (less than 1% of the basin area), the contribution of glacier melting to river discharge in summer and autumn is quite low.In summer (from June to August), the precipitation is 281.6 mm, which is larger than the total of the actual ET (180.1 mm) and runoff (78.3 mm).This result implies that precipitation recharges the soil water and groundwater in summer.In autumn (from September to November), the precipitation is 105.4 mm for the entire basin, and it is lower than the total of the actual ET (58.5 mm) and runoff (62.5 mm).This result implies that the runoff is generated from the precipitation and also from the soil water and groundwater storage.

Spatial Pattern of Water Balance and Relation to Vegetation
To better understand the ecohydrological pattern in the upper Heihe basin, the spatial distributions of the water balance components are further analyzed, including the precipitation, evapotranspiration (ET), runoff and soil moisture in the top one meter during the vegetation growing season from May to October, based on a 30-year model simulation during 1981-2010.Figure 6 shows that annual precipitation and runoff have a similar spatial pattern, annual actual evapotranspiration has a similar pattern with the soil moisture of the top layer in the growing season.The mean soil moisture values are relatively high (0.22-0.41) because of the selected wet season and possible uncertainties of soil water parameters.As shown in Figure 6a, annual precipitation over the study catchment ranges from 220 mm to 630 mm, and the East Tributary has the highest annual precipitation.Comparing Figure 6b with Figure 2, the areas with relatively high ET correspond to the two major vegetation types, namely, steppe and alpine meadow.Figure 6b illustrates that ET has a similar pattern in the East and West Tributaries corresponding with the distribution of alpine meadow vegetation.However, the steep and narrow valley area along the main stream (i.e., downstream of the junction of East and West Tributaries) has relatively low ET and low soil moisture due to the lower precipitation, and this is also related to the steep hillslope and sparse vegetation.In general, the spatial distribution of runoff is mainly controlled by the precipitation and also affected by the topography and vegetation.
The area percentage of each vegetation type is calculated at a 200-m elevation interval in the basin.Accordingly, the water balance components (i.e., annual precipitation, ET and runoff) are also calculated in the same elevation intervals.Figure 7 shows the distributions of vegetation types, annual precipitation, ET and runoff along with the elevation.The vegetation areas of the four major vegetation types (namely, steppe, shrub, coniferous forest and alpine meadow) increase with altitude in the elevation range of 1700-3000 m.Vegetation coverage increases with elevation in the order of steppe, shrub, coniferous forest and alpine meadow.Desert dominates the lowest elevation area, and alpine sparse vegetation dominates the highest elevation area.This spatial vegetation pattern is closely related to the changes in precipitation and temperature along the elevation.Precipitation and runoff have a similar pattern, both increasing with elevation.The actual evapotranspiration shares a similar pattern with the four major vegetation types (i.e., steppe, shrub, coniferous forest and alpine meadow) along the elevation.Actual evapotranspiration increases with elevation when the elevation is lower than 3000 m.The highest actual evapotranspiration is at the elevations of 3000-3600 m, where shrub and alpine meadow are the dominant vegetation types.Previous research has found that during the growing season, vegetation cover is the densest in the   Precipitation and runoff have a similar pattern, both increasing with elevation.The actual evapotranspiration shares a similar pattern with the four major vegetation types (i.e., steppe, shrub, coniferous forest and alpine meadow) along the elevation.Actual evapotranspiration increases with elevation when the elevation is lower than 3000 m.The highest actual evapotranspiration is at the elevations of 3000-3600 m, where shrub and alpine meadow are the dominant vegetation types.Previous research has found that during the growing season, vegetation cover is the densest in the Precipitation and runoff have a similar pattern, both increasing with elevation.The actual evapotranspiration shares a similar pattern with the four major vegetation types (i.e., steppe, shrub, coniferous forest and alpine meadow) along the elevation.Actual evapotranspiration increases with elevation when the elevation is lower than 3000 m.The highest actual evapotranspiration is at the elevations of 3000-3600 m, where shrub and alpine meadow are the dominant vegetation types.Previous research has found that during the growing season, vegetation cover is the densest in the elevation range of 3200-3400 m in the Qilian Mountains [44].This implies that vegetation dynamics are most intensive in the elevation range of 3200-3400 m.When elevation is higher than 3400 m, actual evapotranspiration gradually decreases and vegetation dynamics gradually weaken due to the decreasing temperature.The region with elevation above 4200 m has the highest runoff depth due to high precipitation and lower evapotranspiration.Glacier melting is also a reason for the high runoff depth.Regarding the annual water balance characteristics, the regions with elevations lower than 3200 m are water limited.In these regions, actual evapotranspiration and vegetation growth are controlled by water availability (precipitation).The regions with elevation higher than 3400 m are temperature/energy limited.In these regions, actual evapotranspiration and vegetation growth are mainly controlled by temperature, and the actual evapotranspiration decreases with increasing elevation.Vegetation grows best and the actual evapotranspiration has the highest value in the region with elevations ranging from 3200 to 3400 m.The above spatial patterns of vegetation and water balance components confirm that climate variability has a significant effect on ecohydrological patterns in the high mountainous region.
The water balance of the three different elevation intervals of the East Tributary and West Tributary are calculated and shown in Table 4.The water balance (actual ET and runoff) differs among vegetation types in each elevation interval where the precipitation is close.Both of the runoff depth and runoff ratio for the coniferous forest is less than other vegetation types within the same elevation interval.This implies that vegetation type enhanced the differences in water balance in addition to the climate variability.The water balance of the entire catchment for each vegetation type is also calculated and shown in Table 5, which is the result of the combined effects of climate and vegetation.Because the annual precipitation increases with elevation and different vegetation types grow at different elevations, the annual precipitation of alpine meadow, alpine sparse vegetation and shrub are 488.5 mm, 547.3 mm and 495.9 mm, respectively, which are higher than those of coniferous forest (402.1 mm) and steppe (396.7 mm).The annual average actual evapotranspiration (ET) of coniferous forest, shrub, steppe and alpine meadow ranges from 331.5 mm to 355.0 mm, whereas the actual ET of alpine sparse vegetation is relatively lower (237.2mm).Table 5 shows that the annual runoff depth varies for different vegetation types.Alpine meadow and alpine sparse vegetation have higher annual runoff (147.8 mm and 310.1 mm, respectively) than forest and steppe (70.5 mm and 65.2 mm, respectively).The runoff ratios of the four major vegetation types, namely, steppe, shrub, coniferous forest and alpine meadow, range from 0.16 to 0.30, whereas alpine sparse vegetation has higher runoff ratio of close to 0.5 due to the high altitude.The water yield per unit area from different vegetation type is in order of alpine sparse vegetation, alpine meadow, shrub, coniferous forest and steppe.The top three largest vegetation areas are covered by alpine meadow (with an area of 4549 km 2 ), alpine sparse vegetation (with an area of 2009 km 2 ) and shrub (with an area of 1652 km 2 ), which are also located in relatively high-elevation regions.The total runoff amount produced by the top three largest vegetation areas are 6.72 ˆ10 8 m 3 /a (alpine meadow), 6.23 ˆ10 8 m 3 /a (alpine sparse vegetation) and 2.33 ˆ10 8 m 3 /a (shrub) (see Table 5).The runoff amount produced by forests (with area of 561 km 2 ) is 0.40 ˆ10 8 m 3 /a.The coniferous forest has a small area (561 km 2 ) and produces a small amount of runoff (0.40 ˆ10 8 m 3 /a).
Figure 8 shows the decadal changes in water balance components.Annual precipitation, actual ET and runoff show similar patterns along the elevation.However, decadal variability of the water balance is also observed.Compared with the 1980s, precipitation decreased and ET increased in the 1990s, which caused a reduction in runoff.However, in the 2000s, precipitation increased in the high-altitude region with elevations above 3600 m, which led to a significant increase in runoff compared to the 1990s.The water balance of the entire catchment for each vegetation type is also calculated and shown in Table 5, which is the result of the combined effects of climate and vegetation.Because the annual precipitation increases with elevation and different vegetation types grow at different elevations, the annual precipitation of alpine meadow, alpine sparse vegetation and shrub are 488.5 mm, 547.3 mm and 495.9 mm, respectively, which are higher than those of coniferous forest (402.1 mm) and steppe (396.7 mm).The annual average actual evapotranspiration (ET) of coniferous forest, shrub, steppe and alpine meadow ranges from 331.5 mm to 355.0 mm, whereas the actual ET of alpine sparse vegetation is relatively lower (237.2mm).Table 5 shows that the annual runoff depth varies for different vegetation types.Alpine meadow and alpine sparse vegetation have higher annual runoff (147.8 mm and 310.1 mm, respectively) than forest and steppe (70.5 mm and 65.2 mm, respectively).The runoff ratios of the four major vegetation types, namely, steppe, shrub, coniferous forest and alpine meadow, range from 0.16 to 0.30, whereas alpine sparse vegetation has higher runoff ratio of close to 0.5 due to the high altitude.The water yield per unit area from different vegetation type is in order of alpine sparse vegetation, alpine meadow, shrub, coniferous forest and steppe.The top three largest vegetation areas are covered by alpine meadow (with an area of 4549 km 2 ), alpine sparse vegetation (with an area of 2009 km 2 ) and shrub (with an area of 1652 km 2 ), which are also located in relatively high-elevation regions.The total runoff amount produced by the top three largest vegetation areas are 6.72 × 10 8 m 3 /a (alpine meadow), 6.23 × 10 8 m 3 /a (alpine sparse vegetation) and 2.33 × 10 8 m 3 /a (shrub) (see Table 5).The runoff amount produced by forests (with area of 561 km 2 ) is 0.40 × 10 8 m 3 /a.The coniferous forest has a small area (561 km 2 ) and produces a small amount of runoff (0.40 × 10 8 m 3 /a).
Figure 8 shows the decadal changes in water balance components.Annual precipitation, actual ET and runoff show similar patterns along the elevation.However, decadal variability of the water balance is also observed.Compared with the 1980s, precipitation decreased and ET increased in the 1990s, which caused a reduction in runoff.However, in the 2000s, precipitation increased in the highaltitude region with elevations above 3600 m, which led to a significant increase in runoff compared to the 1990s.In this study, the model simulation shows that the precipitation recharges soil water and groundwater in the summer.This result is consistent with a previous study by Yang et al. [72] on the Heihe River.Using the hydrochemistry approach, it is found that the precipitation contributes only slightly to the total surface runoff, instead mainly recharging the sub-surface soil and groundwater.The interflow and groundwater flow dominate the total runoff.Additionally, Wang et al. [73] inferred that as the air temperature has risen in the Heihe basin during the last 50 years, the active layer of frozen soil has increased in thickness.This may lead to the increase of both soil water storage and interflow and thus change the contribution of each runoff component in the future.
Our study shows that forests contribute only a small amount of the water yield.In a previous study, He et al. [3] analyzed the water balance of a small experimental catchment (2.91 km 2 ) in the upper Heihe basin and found that the forests contribute only 3.5% to the total runoff, which is similar to this study.This result is also supported by model simulations.For example, using the Topography Driven Flux Exchange (FLEX-Topo) model, Gao et al. [74] report that forest hillslope generates only a small amount of runoff in the upper Heihe basin.Qin et al. [75] analyze the water balance components of different landscapes in the upper Heihe basin using the VIC model and find that glaciers contribute 3.57% of the total runoff and that the contribution of forests is also quite small (0.49%).This result also implies that the barren regions contribute most of the total runoff (52.46%) and steppe contributes 34.15% to the total runoff.This is in accordance with our results, as Qin et al. [75] consider the steppe and alpine meadow as the same vegetation type and consider regions in high elevation as barren.
Our study shows that the runoff ratio of coniferous forest has the lowest value (0.18) comparing with other vegetation types except the desert in the upper Heihe River basin.This is similar to the findings by Yaseef et al. [76], which shows that when precipitation larger than 300 mm ET accounts 85% of precipitation for Aleppo pine forest in Southern Israel.Wang et al. [77] analyzed the annual water balance of shrub at a station in the Inner Mongolian Highland Region with elevation of 1300 m and found that ET/P ratio is about 94%, which is higher than the results of the present study because of the differences of elevation and air temperature in the two study areas.Wang et al. [78] analyzed the water balance of different vegetation in a small watershed in the Liupan Mountains, Northwest China, based on field measurements.They found that the evaporation rate (ratio of ET to precipitation) was about 60% for grassland, 93% for shrubs, and >95% for forest.This also shows that the forest has the least water yield compared with other vegetation types in the Liupan Mountains.

Conclusions
A geomorphology-based ecohydrological model (GBEHM) is developed in the upper Heihe basin, and this model is validated using available observations, including soil moisture, streamflow discharge, and actual evapotranspiration estimated from remote sensing.The catchment water balance characteristics and their spatial-temporal variability are analyzed based on the ecohydrological simulation.The following conclusions can be drawn from the results of this study: (1) At the basin scale, the model provides a good simulation of streamflow discharge in the two tributaries and the entire catchment of the study area.It also captures the spatial pattern of soil moisture appropriately.In addition, the simulated actual evapotranspiration and remote sensing-based estimation have close long-term average values and similar spatial patterns over the entire study catchment.The GBEHM may be useful for ecohydrological simulation and prediction in cold high-altitude regions.(2) Analysis of the water balance characteristics shows that water balance characteristics are closely related to the altitude and vegetation patterns in the study catchment.Regarding the annual water balance characteristics, the low-altitude regions with elevations below 3200 m are water limited.
The actual annual evapotranspiration and vegetation distribution and growth are controlled by water availability (precipitation).Seasonal analysis indicates that river runoffs are mainly in summer and autumn, and runoff in spring is generated from precipitation and snow melt.(3) In the upper Heihe basin, the precipitation and runoff share a similar pattern, increasing with elevation.Actual evapotranspiration has a similar pattern with the four major vegetation types (i.e., steppe, shrub, coniferous forest and alpine meadow) along the elevation.The highest actual evapotranspiration is at the elevations of 3000-3600 m, where shrub and alpine meadow are the two dominant vegetation types.Precipitation controls the spatial pattern of annual runoff and determines the spatial pattern of vegetation together with the air temperature.Climate variability in the high mountainous region has a significant effect on ecohydrological patterns.( 4) At the same time, vegetation type enhanced the differences in annual runoff and actual evapotranspiration.In the same elevation interval with similar precipitation, differences in the runoff depth (and the actual evapotranspiration) were caused mainly by the vegetation types.
For the whole study area, the water yield per unit area from different vegetation types is in order of alpine sparse vegetation, alpine meadow, shrub, coniferous forest and steppe.The three major vegetation types, namely, alpine meadow (with an area of 4549 km 2 ), alpine sparse vegetation (with an area of 2009 km 2 ) and shrub (with an area of 1652 km 2 ), located in relatively higher elevation contribute most of the river runoff.
Several limitations remain in the current study.The GBEHM simulates the vegetation dynamics with the known leaf area index and other vegetation parameters.Further improvement of the model should include carbon partitioning to simulate the vegetation growth.The uncertainty of the model parameters should be assessed to apply this model to other catchments.

Figure 1 .
Figure 1.Location and topography of the upper Heihe basin.

Figure 1 .
Figure 1.Location and topography of the upper Heihe basin.

Figure 2 .
Figure 2. Vegetation map of the upper Heihe basin.

Figure 3 .
Figure 3.Comparison of the monthly observed and simulated river discharge at the three hydrological stations: (a) Yingluoxia; (b) Zhamashike; and (c) Qilian, in the 2001-2012 period.

Figure 3 .
Figure 3.Comparison of the monthly observed and simulated river discharge at the three hydrological stations: (a) Yingluoxia; (b) Zhamashike; and (c) Qilian, in the 2001-2012 period.

Figure 4 .
Figure 4. Spatial distributions of monthly average soil moisture (m 3 •m −3 ) in the East Tributary: (a) observed soil moisture at 4 cm from the surface and (b) simulated soil moisture of the top 5-cm layer.

Figure 4 . 21 Figure 5 .
Figure 4. Spatial distributions of monthly average soil moisture (m 3 ¨m´3 ) in the East Tributary: (a) observed soil moisture at 4 cm from the surface and (b) simulated soil moisture of the top 5-cm layer.Forests 2016, 7, 10 11/21

Figure 5 .
Figure 5.Comparison of the annual average evapotranspiration between the (a) GBEHM simulation and (b) remote sensing-based estimation during the 2001-2012 period [59].

Figure 6 .
Figure 6.Spatial distributions of the mean values of (a) annual precipitation; (b) annual evapotranspiration; (c) annual runoff and (d) growing season soil moisture of the top 1 m during 1981-2010.

Figure 8 .
Figure 8. Changes in annual average precipitation, evapotranspiration and runoff along with elevation during the 1980s, 1990s and 2000s.

Figure 8 .
Figure 8. Changes in annual average precipitation, evapotranspiration and runoff along with elevation during the 1980s, 1990s and 2000s.

Table 1 .
Major parameters related to vegetation and soil types in the geomorphology-based ecohydrological model (GBEHM).

Table 2 .
Water balance for the East and West Tributaries and the entire catchment during the 2001-2012 period.

Table 2 .
Water balance for the East and West Tributaries and the entire catchment during the 2001-2012 period.

Table 3 .
Seasonal water balance for the East and West Tributaries and the entire catchment during the 2001-2012 period.

Table 4 .
Water balance in different elevation intervals of the East Tributary and West Tributary during the 2001-2012 period.

Table 5 .
Water balance of the entire catchment for each vegetation type during the 2001-2012 period.