Vegetation response to upstream water yield in the Heihe River by time series analysis of MODIS data

L. Jia, H. Shang, and M. Menenti State Key Laboratory of Remote Sensing Science, Jointly Sponsored by the Institute of Remote Sensing Applications of Chinese Academy of Sciences and Beijing Normal University, Beijing 100101, China Alterra, Wageningen University and Research Centre, Wageningen, The Netherlands Wageningen University, Centre for Geo-Information, Wageningen, The Netherlands Department of Earth Observation, Delft University of Technology, Delft, The Netherlands


Introduction
The Heihe River is the second largest inland river in the northwest (NW) China, it originates from the Qilian Mountains in Qinghai province, runs through the Hexi Corridor of Gansu Province and flows into the western Inner Mongolian Plateau (Fig. 1).Liquid and solid precipitation is abundant in the high elevation, upper reach of the river, which produce sufficient runoff going to the down streams.The mid and lower reaches of the Heihe River, located in an arid and semi-arid climate environment, are characterized by very limited precipitation.The economy, agriculture and human living conditions in the Heihe River Basin are therefore highly dependent on the streamflow of the Heihe River.Streamflow from the upper catchment of the basin becomes the major water resource besides the ground water.Ground water, replenished by local precipitation and water flow from the upstream catchment, also plays an important role in the ecohydrological cycle of the low reach.However, in the last thirty years, due to the large increase of population and the fast development of economy, the water competition between human activities and nature became more and more serious, particularly in the middle and lower reaches.The development of modern irrigation schemes in the middle reach of the basin is taking up an increasing share of fresh water resources, endangering the oasis and traditional irrigation systems in the lower reach, the Ejina oasis (Fig. 1).The natural vegetation in this area has significantly decreased since 1990s and the desertification became more and more serious in the lower reaches of the Heihe River.A better understanding of the relationship between vegetation growth and the water availability could help managing the limited water resources more efficiently and protect the vulnerable biosphere along the Heihe River.
Satellite remote sensing data has been successfully applied to monitor and study vegetation dynamic and land cover changes (Bewket, 2002;Justice and Hiernaux, 1986;Zomer et al., 2001).The Normalized Difference Vegetation Index (NDVI) assumes that the difference between the near-infrared and red reflectance divided by the sum of both is a quantitative measure of photosynthetic activity.The NDVI is a quantita- tive description of vegetation conditions directly obtained from satellite remote sensing data and is widely used in many different application fields, for example, drought assessment (e.g.Liu and Negron Juarez 2001), agricultural productivity assessment (e.g.Groten and Ocatre 2002;Maselli and Rembold 2001), climate linkages (e.g.Roerink et al., 2003;Wen et al., 2004) and so on.Recently, with the growing collections of satellite remote sensing data, the analysis of NDVI time series data has paved the way towards the study of phenological variations during a relatively long time period.
In order to understand the evolution of vegetation in response to water availability and the causes of phenological variations, a high temporal resolution is preferable in the NDVI time series.Optical remote sensing measurements are influenced by the atmosphere, clouds, bi-directional effects, the solar elevation angle and so on, thus the time series of NDVI data contain erroneous observations which need to be identified and removed.This yields gaps either in space or in time.Cloud contamination is the major factor affecting NDVI since lower NDVI are observed under hazy conditions and very low or negative values are observed when cloudy.Such cloud-contaminated observations are not easily identified and removed from the original data.Cloudiness not only reduces the reliability of NDVI data, but also makes it unfeasible to obtain a continuous NDVI time series for a long time over a large area.In addition, inadequate atmospheric correction might introduce noise, which would lead to errors in data analysis and interpretation.To reduce the occurrence of cloud-contaminated observations, most NDVI data sets are Maximum Value Composite (MVC) products (Holben, 1986), such as the Pathfinder AVHRR dataset and the SPOT VGT dataset.Such data products are generated by selecting and retaining the maximum NDVI value observed in each pixel during a certain period.This method can significantly improve the quality of the NDVI time series, but it cannot totally remove the noise due to cloud-contamination.
To reduce the noise and generate a cloud-free NDVI time series for data analysis, a number of methods have been developed and have been successfully applied in some cases.Some examples include the Best Index Slope Extraction algorithm (BISE) (Viovy et al., 1992;Lovell and Graetz, 2001) for analysis of daily NDVI data, the Asymmetric Gaussian function (Jonsson and Eklundh, 2002), a Savitzky-Golay filtering algorithm (Chen et al., 2004), the Mean-Value Iteration filter (MVI) (Ma and Veroustraete, 2006a), and the weighted least-squares linear regression to the temporal NDVI signal (Swets et al., 1999), among others.Though these algorithm have been successfully used in different applications, limitations are also existing.For instance, the length of the sliding window and threshold value in BISE method may need to be adapted to different phenological stages and for different plant species.The weighted least-squares linear regression and the Asymmetric Gaussian function could not capture the complex phenology of land cover characterized by two or more growth cycles in one year.The reliability of the Savitzky-Golay filtering algorithm strongly depends on the assumption that the envelope of the original data gives the best description of vegetation growth, while this concept may lead to overestimate of the NDVI values, which are in most cases generated using the MVC method.When two or more adjacent cloud contaminated NDVI values occur, many iterations with the MVI method are required to estimate a reliable NDVI value.
Besides the methods mentioned above, the Discrete Fourier Transform was also used in NDVI time series analysis.The signal for each individual pixel can be decomposed into its Fourier components, such as series of harmonic sine or cosine waves, from which a cloud-free signal can be reconstructed (Menenti et al., 1993;Menenti et al., 1995;Verhoef, 1996;Azzali and Menenti, 2000).In the Fourier transform of NDVI time series, the term harmonic implies that the time series is periodic.The assumption that the natural variation can be reproduced by a periodic model was inspired by the fact that time series of yearly and daily datasets of remotely sensed images often clearly displayed the influence of the diurnal and yearly courses of the sun on the earth's surface.The actual signals were not purely periodic, but the periodic model could be successfully applied to obtain insight in the dynamics of land surface processes.So it was widely used to study the linkages between climate and vegetation dynamic (e.g.Roerink et al., 2003;Wen et al., 2004;Immerzeel et al., 2005;Julien et al., 2006).Different algorithms were developed, i.  , 1993) and Harmonic ANalysis of Time Series (HANTS) (Verhoef, 1996).The HANTS algorithm was developed to deal with time series of irregularly spaced observations and to identify and remove cloud-contaminated observations.The Fourier analysis problem was formulated as a curve fitting problem, so that it was possible to reject certain observations by attaching a weight of zero to them.The parameters describing the curve were determined from the weighted least square fitness.Once the model is determined using the reliable (cloud-free) observations only, it can be used to gaps generated by removing the cloud-contaminated observations.The objective of this research was to study the vegetation dynamics in the Ejina Oasis, located in the terminal zone of the lower reach of the Heihe Basin in response to water availability, including both local rainfall and water flow from the upper and mid reaches, with particular attention to the significant inter-annual and intra-annual variability.Since the Ejina Oasis is located in the inland delta of the Hei He, different hydrological processes contribute to determine vegetation phenology: local but scarce rainfall, water flow through the channel system in the delta and water flow through the shallow unconfined aquifer.As regards the latter, further complexity is added by the role of capillary rise and root water uptake, which both depend on soil type and on the depth of the shallow groundwater.We used one Fourier-Transformation method -Harmonic Analysis of Time Series (HANTS) algorithm -to reconstruct time series of cloud-free NDVI data and analyze the vegetation responses to the water availability in the lower reaches of the Heihe River.

Study area
The study area is the arid and semi-arid land in Ejina Oasis, located in the downstream of Heihe River Basin in the northwest China (Fig. 1), across Qinghai, Gansu  (Wang and Cheng, 1999), mean annual pan evaporation more than 3500 mm (Jin et al., 2008), and annual mean sunshine time duration of 3000 to 4000 h (Ma and Veroustraete, 2006b).Hence, the natural vegetations are fragmented and easily disturbed in this area.Fig. 2 shows the land cover in the Ejina Oasis according to the data set provided by the Environmental and Ecological Science Data Center for West China, National Natural Science Foundation of China (http://westdc.westgis.ac.cn) in 2004.Shrubland and grassland are the major land cover types in Ejina Oasis.There is some dry cropland near the Heihe River.The grassland and the cropland are the main economic resource in the Ejina Oasis.The forest land consists of Poplar P. Euphratica.The Heihe River is the only water resource in the Ejina Oasis, which flows from the Qilian mountain, through the Hexi Corridor and flows into Inner Mongolia.The Yingluoxia gauge station and the Zhengyixia gauge station (Fig. 1) are located at the upper boundary of the mid and lower reach of the Heihe River, respectively.The river divides into an Eastern and Western branches past the Langxinshan station.The water flow reaching the lower reach is strongly dependent on the water consumption in the mid reach, where the irrigation system is well developed.After the construction of water reservoirs and water diversion infrastructure, water supply to the Ejina Oasis has been significantly reduced.Regulation of river flow through the operation of reservoirs led to the degeneration of vegetation and desertification in the Ejina Oasis (Qi and Luo, 2006).

Satellite data
In this research, we used the Terra-MODIS VI dataset (MOD13A2), which covers the period from 2000 to 2008, with a spatial resolution of 1 km.The NDVI data was pro- duced from the Moderate Resolution Imaging Spectroradiometer, successor to the AVHRR instrument, but with different red and NIR channels, i.e. 620-670 nm and 841-876 nm.The data product is generated using the MVC algorithm over each subsequent 16-days period, starting from the first day of the year.

Meteorological and hydrological data
Besides the remote sensing data sets, we also used meteorological and hydrological data from the Heihe River Basin.There are five meteorology stations in the Heihe River region (Fig. 1), which are Dingxin, Gaotai, and Zhangye in Gansu Province, Ejina Banner in Inner Mongolian Province, and Qilian in Qinghai Province.The Ejina Banner station is located in Ejina Oasis and was considered to represent the precipitation conditions in the study area, while Dingxin is located at the start of the down stream (just after the Zhengyixia station along the river) (Fig. 1).Data from these two stations will be used in this study.The monthly precipitation was calculated from daily measurements at the two meteorological stations from 1952 to 2007.
The hydrological data was collected from two hydrological stations located at the start and at the end of the mid reach, just upstream of the Ejina Oasis (Fig. 1).The monthly streamflow at each station was calculated from daily measurements from 1978 to 2004.

HANTS Algorithm
The HANTS algorithm (Menenti et al., 1993;Verhoef et al., 1996;Roerink et al., 2000) is based on the discrete Fourier Transform.y is given, and y i indicates the i -th sample, then the series can be described by means of a Fourier series: where ω j is the frequency of the j th -harmonic term of the Fourier series, t i is the time at which the i -th sample was taken and M is the number of the frequencies of the Fourier series (M < = N), a j and ϕ j are the amplitude and the phase of the j th -harmonic term, respectively.As the 'zero' frequency has no phase and the cosine becomes one, the amplitude associated with the zero frequency, a 0 , is equal to the average signal of the N samples of variable y.
In the HANTS algorithm the most important control parameter is the frequency ω j which is defined according to the most significant frequencies expected to be present in the time series of observations.After setting the number of frequencies M and frequencies ω j , the unknown parameters of the Fourier series are the amplitudes a j and the phases ϕ j values, which are determined by fitting the time series of observations.HANTS handles the Fourier analysis as a curve fitting problem.To evaluate the effectiveness of any observation in the time series, the HANTS was designed by attaching different weights to different observations according to the influence of cloudcontaminated observations.In the current version of HANTS, only two weights options are made available -either "one" (corresponding to "good" data) or "zero" (corresponding to "bad" data).To remove the negative influence of cloud-contamination on the NDVI value, the HANTS performs curve fitting iteratively: in the first step, the least squares curve fitting is performed using all data.In the second step, the observations which are below the fitting curve and whose distance is larger than a pre-defined distance are removed by attaching a weight of zero.The remaining data are used to compute the least squares curve fitting again, then outliers are identified as in the first step.This iterative procedure is repeated until all the remaining observations are within a pre-defined distance of the fitting curve.4185

Implementation of HANTS and parameters determining performance
To obtain a reliable fitting curve, seven essential parameters have to be specified to run the HANTS algorithm: -Valid data range -the range of the input variables values.The values out of this range are assigned a weight of zero.
-Period -the temporal length of each term in the Fourier series; -Number of Frequencies (NOF) -the number of harmonic terms over one period.
The complex fluctuation of the time series data might need more frequencies to get an accurate fit.The zero frequency is automatically considered in the algorithm as the average value of the valid observations in the time series.
-Direction of outliers -is used to indicate the direction of unreasonable values with respect to the current curve.In the HANTS algorithm, these unreasonable values are called "outliers".In time series of NDVI observations cloud-contaminated pixels usually have low NDVI values, the direction of 'outliers' is set to "low" which means that a NDVI value in the time series smaller than the predefined threshold will be eliminated by giving a weight "zero" to it.-Fit Error Tolerance (FET) -is the distance which identifies observations in the times series farther than this distance as "outliers".If all the remaining observations are within a distance of FET from the fitting curve, the iterations will stop.
-Degree of Over Determinedness (DOD) -the minimum number of extra data points, which have to be used in the ultimate fit.The number of valid observations must always be greater than or equal to the number of parameters that describe the curve (defined by 2× NOF−1).In order to get a more reliable fit more data points than the necessary minimum should be included as specified by the DOD value.There is no straightforward rule in the determination of these essential parameters, some practice is usually needed before getting optimal values.More application cases can be found in literature (e.g.Azzali and Menenti, 2000;Roerink et al., 2000;Wen et al., 2004 among others).

Temporal-similarity-statistical (TSS) method to improve HANTS performance
Besides a proper setting of essential parameters described above, the reliability of the reconstructed NDVI time series also depends on the quality of the original time series of NDVI data.The quality of the observation time series is associated to, for instance, probability of missing observations, number of gaps, size of continuous gaps in the time series, location of the gaps, and retrieval quality of land surface parameters.Large probability of missing data, especially the large size of continuous gaps may result in an unreasonable reconstructed curve, which limits the application of the HANTS algorithm.
A two-step method was developed in this study to deal with long gaps, which is called temporal-similarity-statistical (TSS) method.The method is based on historical data and finds in the available data record the initial values to fill extended gaps.
The algorithm compares the current year time series with the available data record to construct segments of the time series to fill extended gaps.
For each pixel, a temporal coordinate system for the time series of NDVI in one year growth cycle including all sample composites (for instance 23 composites for MODIS VI products with 16-day interval) is built first as shown in Fig. 3.In Fig. 3, The horizontal axis stands for the composite coordinate in a given year, while the vertical axis for year coordinate of the same composite.
where CV i j is the composite variance between the i th year and the j th year, a ki is the original NDVI of the kth composite in i th year, A kj is the reconstruct NDVI of the kth composite in the j th year, N is the number of the reference samples used to calculate the time series of reconstruct data in i th year.
The year having the smallest CV value will be used as the reference year, the composite in the reference year time series corresponding to the missing data period of the target year will be taken as the initial estimates of these missing observations.This procedure is only applied when the gap is so big that will affect the reconstructed curve when using the standard HANTS algorithm.The next step is therefore to evaluate the size of the large gaps in the time series that have significant influence on the reconstruction results of NDVI time series.The threshold gap size is found to be correlated to the length of term of the highest frequency of the Fourier transform component.
The days of missing data should be less than half of the length of the highest frequency term of the Fourier series.

Intensity-Hue-Saturation (IHS) transform
To analyze the vegetation growth condition especially the phase information, the IHS (Intensity, Hue, Saturation) transform was developed in the HANTS algorithm.In the HIS transform, for each pixel, the mean value is assigned to the intensity component, the phase of to the color hue, and the amplitude to the color saturation.In HANTS algorithm, the transform between HIS and RGB (red/green/blue) is carried out through the following formulas (Verhoef et al., 1996;Julien et al., 2006): Here r, g and b are the color signals in red, green and blue, M is the mean signal, A is the amplitude, P is the phase in degrees.The mean signal is scaled by the minimum value C and a saturation value S. The amplitude is scaled between zero and a maximum value A max When the period follows one year cycle, the phase goes from 0 • and 360 • (say 1 • corresponding to about 1 day).
In this study, the amplitude of the "zero frequency" of the Fourier series, e.g. the annual mean NDVI, was assigned to the intensity component.The amplitude of the 12-month term is found to be much larger in the Ejina area than those of the other three terms, the phase of the 12-month term can be taken as the timing of the peak green-up and correlated to the color hue.The amplitude of the 12-month harmonic term is assigned to the color saturation.The order of the colors is therefore bluecyan-green-yellow-red-magenta-blue, corresponding to the phase timing in the month of January-March-May-July-September-November-January, respectively.

Anomaly calculations
The anomaly of annual precipitation and annual rive streamflow was defined by where k is the year number, ∆V k is the anomaly of annual precipitation or river streamflow in the year k, V k is the annual precipitation or river streamflow in the year k, N i =1 V i is the sum of annual precipitation or river streamflow between the N years.
The monthly anomaly of streamflow was calculated as where j is the month number and k is the year number, ∆Q k,j is the month j anomaly in year k, Q k,j is the monthly total streamflow in month j and year k, and N is the mean monthly streamflow for month j averaged over N years (between 1980 and 2004) for the same month.
To calculate the yearly anomaly in NDVI, V k in Eq. ( 4) is replaced by annual mean value of NDVI.

Precipitation
In the period from 1980 to 2007 the annual precipitation did not exceed 100 mm at Dingxin and 80 mm at Ejina Banner (Fig. 4a).The mean annual precipitation was 54 mm at Dingxin and 32 mm at Ejina Banner.Anomalies calculated by Eq. ( 3) for each year's annual precipitation are presented in Fig. 4b.The annual precipitation in Ejin Banner showed a periodic oscillations in the last 28 years.From 1980 to 1989, the annual precipitation at Dingxin and Ejina Banner was lower than the long term average in most of the years.From 1990 to 1999, more annual precipitation than the previous ten years was observed both at Dingxin and at Ejina Banner, in particularly in 1995 and 1999.Since 2000, the annual precipitation in Dingxin and Ejina Banner decreased significantly again, lower than the long-term average.
Taking into account that mean annual potential evaporation in the same period was 3500 mm (Jin et al., 2008), we may conclude that the contribution of precipitation to water supply of vegetation in the Ejina Oasis is negligible.

Streamflow
The annual streamflow at Yingluoxia station, located at the start of the middle stream area of the Heihe River region, stayed at a relative high level with the 23-years average of 1.66×10 9 m 3 (streamflow data in 2001 and 2002 at Yingluoxia station were not available) (Fig. 5).However, the annual streamflow at Zhengyixa station, located at the end of the middle stream area, showed a significantly decrease since 1990.Before 1990, the average of the annual streamflow at Zhengyixia was 1.11×10 9 m 3 , while it decreased to 0.75×10 9 m 3 from 1990 to 2001 due to a significant increase in water demand in the middle reach areas of the Heihe River.The streamflow from the Heihe River supplied most of water consumed by the ecosystem in Ejina Oasis.The large decrease of streamflow from middle stream area caused many environmental problems in the lower reach of the Heihe River Basin , such as oasis desertification and vegetation degeneration (Qi and Luo, 2006).The anomaly of streamflow for each year was calculated similar to the annual precipitation anomalies by Eq. ( 2).The anomalies at the two stations showed a synchronous increase or decrease in most years, but with different values, which indicated that the streamflow at Zhengyixia station highly depended on the streamflow at Yingluoxia station and the water consumption in the middle stream area (Fig. 5b).
The monthly streamflow at Zhengyixia station in Fig. 6a showed that the discharge water from January to March was stable from 1998 to 2004, higher than 7×10 7 m 3 .A sharp decrease in the streamflow was observed from April, reached its minimum in May and started to increase again in July.From May to June there was a cut-off of water supply to the Ejina Oasis from the Heihe River, this is mainly because irrigation started in the agricultural area of the middle reach.The river had no continuous flow to the Ejina Oasis for more than 100 days per year since 1950 (Zhang et al., 2005).
From July, the water supply from Zhengyixia station down to the Ejina Oasis increased again and reached the peak in the next three months.The streamflow decreased from October and increased in December again.Large change of monthly streamflow mainly occurred between July and October as shown by the anomaly of the monthly streamflow at Zhengyixia station in Fig. 6b.In 1999In , 2000In , and 2001, the anomaly of the streamflow in most months between July and October showed large negative values, implying that the vegetations in the Ejina Oasis might suffer from lack of water supply in these three years.

Annual mean area-averaged NDVI
Vegetation in the Ejina Oasis depends critically on river flow released in the upper and mid reaches, past allocation and use for irrigation.The operation of reservoirs and of the irrigation infrastructure in the mid reach leads to a large intra-annual variability in river flow released to the lower reaches (Fig. 6), characterized by sudden, short-lived increases in river flow.Large annual variability was characterized by continuous few years negative annual anomaly (Fig. 5).In this section, the annual vegetation growth conditions were analyzed by linking to the inter-annual and intra-annual fluctuations of the water yield from the upper stream of the Ejina oasis.The amplitude of "zero frequency" of the Fourier components of the NDVI time series covered shrubs and sparsely covered grassland, referred to as "oasis" in the following analysis.Desert and Gobi pixels were not considered.The area-average of annual mean NDVI of Ejina Oasis was therefore calculated by averaging the annual mean NDVI over the "oasis" pixels in the downstream area.Taking the Ejina Oasis as a whole, the difference between the drier conditions in the years 2000 through 2002 and wetter conditions in 2003 and 2005 is clearly reflected by the evolution of the annual mean NDVI (Fig. 7).The vegetation condition was the worst in 2001, while it was gradually improved with the best condition appeared in 2004 (Fig. 7).
The timing of the vegetation green-up, i.e. peak NDVI, is one of the indictors of vegetation growth with respect to the influence of climate and water availability.The amplitude of the 12-month harmonic component was found much larger than the amplitudes of the 6-month and the 3-month harmonic components, the 12-month harmonic component therefore reflects the main growth cycle of the vegetation in the Ejina oasis.Hence, the phase value of the 12-month harmonic component of the Fourier Transform gives the timing of the peak green-up of the vegetation.A significant delay in the timing of peak green-up of the vegetation, characterized by a larger phase value of the 12-month harmonic, occurred in 2001, the driest year in the sample observed (Fig. 7).

Influence of water availability and climate on vegetation growth
The time series of river flow and of the area-averaged monthly NDVI (Fig. 8) clearly shows that vegetation green-up trails river flow, but also shows no response of vegetation green-up to the intra-annual fluctuations in river flow.The estimation of the time lag between river flow and vegetation green-up and especially its inter-annual variability is of importance towards a better understanding of hydro-ecological processes in the Ejina Oasis.
The was found between the annual vegetation condition and the previous-year streamflow (R 2 =0.86,Fig. 9a), while no correlation between the annual mean NDVI and the current year total streamflow was found.Similar results was also found by Jin et al. (2008).However, the correlation with the streamflow of the current year's spring was significant (R 2 =0.96,Fig. 9b).Higher annual water yield from the upper stream of the Heihe River in the previous year might have provided sufficient water to recharge the ground water in the down stream of the river, i.e. the Ejina Oasis.Both correlations imply a significant role of shallow groundwater in determining the delayed response of vegetation conditions to water supply.We have analyzed the time lag between streamflow and vegetation green-up by using the phase value of the 12-month harmonic.A plot of the area-averaged phase value versus previous year streamflow reveals two completely different patterns (Fig. 10a): the timing of the peak green-up increased by up to about 25 days (i.e. a 25 • increase in the phase value) over the years (2000,2001,2002) with a slightly negative, but similar anomaly in streamflow, with the largest increase observed in the driest year (i.e.2001).
Under higher than average streamflow, the timing of peak green-up became earlier.Similar results were observed when considering streamflow in the spring of the current year (Fig. 10b).In both cases the 2000, 2003 and 2004 data suggest a positive linear correlation with streamflow.This clearly requires further investigation by extending the concurrent time series of streamflow and NDVI.Here we have evaluated the temperature sum (i.e.accumulated degree-days) over the growing season, which is the main determinant of the duration of phenological stages.The duration of phenological stages is determined by the duration of sustained higher temperature, specifically measured by the sum of temperatures exceeding a species (variety) specific threshold and expressed in degree-days.Here the temperature sum was calculated using a generic threshold of 10 • C. The results show (Fig. 11) that increasingly higher values are obtained under increasingly drier conditions (2001 and 2002) (see also Figs. 4,5 and 7) and that an earlier peak green-up and lower values of the temperature sum were observed under wetter conditions (2003 and 2004).Vegetation development and green-up in 2000 was clearly limited by the rather large temperature sum, the largest observed during the five years considered in this study.
The duration of the vegetative growth stage is particularly sensitive to sustained higher temperatures as measured by the temperature sum.

Spatial variability of vegetation growth conditions
Using IHS color transform described in Sect.3.2, time when maximum NDVI appears of a pixel in different years can be identified.Since the river streamflow data collected at the Zhengyixia gauge station were only available till 2004, the analysis and discussion on the vegetation response to the water availability in Ejina Oasis of the downstream of Heihe River were mainly limited to the period between 2000 and 2005, though IHS maps till 2008 are displayed in Fig. 12.
In an IHS transform color image, light grey to white stands for a nearly constant high NDVI during whole year, while dark refers to a constantly low NDVI, such as desert area.The color hue indicates the appearing time of the maximum NDVI.In Ejina Oasis, the vegetation growth is usually between April and October.Thus if there was vegetation growing in the pixel, the color of the pixel should be in the range of green-yellow-red.The pixels with color out of this range are covered by desert or building area.With the similar hue but different intensive, for example, one pixel with bright yellow, the other with dark yellow, the bright pixel would have a higher annual mean NDVI than the other's, but the time of maximum greenness occurring is close to each other.
In Fig. 12, the area with black or dark blue around the Ejina Oasis is desert or Gobi around Ejina Oasis.Within the Ejina Oasis the peak green-up of most vegetated area appeared among July, August, and September (color from yellow, orange to red) from 2000 to 2008, which covered most areas in Ejina Oasis.The dark green area located between the Ejina Oasis and its surrounding desert area is very sparsely covered area by vegetation.The blue area inside Ejina Oasis is probably the Heihe River or building area.The area with magenta indicates the maximum NDVI appeared in October or November, which could be degraded trees.With the similar hue but different intensity, for example, one pixel with bright yellow, the other with dark yellow, the bright pixel would have a higher annual mean NDVI than the other's, but the time of maximum greenness is close to each other.In 2004 and 2007, most of the Ejina oasis was characterized by a bright yellow color indicating the vegetation condition in these two years were rather good.In the north-west and south-east part of the oasis in 2002, brighter "yellow" color than any other years implies that in these areas vegetation was in the best conditions compared to other years.
The timing of peak green-up did change in different years: In 2000, most pixels in Ejina Oasis had a color of yellow or orange, which implies that the maximum NDVI in most vegetation area occurred in the mid of July or beginning of August.
In 2001, the color of most pixels changed to red, indicating that the timing of maximum greenness was significantly delayed to the end of August or beginning of September.
From 2002 to 2004, number of pixels with yellow and orange were increased gradually, which indicated that the timing of maximum NDVI was earlier than that in 2001, shifted back to August or July.In 2002, "orange" and "red" were dominant colors, while in 2003 the 'orange' was dominant color with significantly increased "yellow" color.In 2004, vegetation condition seems the best expressed by large fraction of yellow pixels, the size of oasis area has also slightly increased when compared to 2000.
Taking the Ejina Oasis as a whole, the sequence of the timing of peak green-up for the 9 years period is estimated as 2000-2007-2004-2008-2003-2005-2002-2006-2001.To evaluate whether the South-North (S-N) propagation of the fluctuations in streamflow through the drainage channels and the shallow aquifer system had any influence on the timing of peak green-up, we calculated the mean phase value in a East-West direction, i.e. over all columns of the 12-month phase image for each row, and plotted such mean value in the S-N direction, i.e. along the rows of the same image (Fig. 13).
A significant gradient appears only in drier years (2001 and 2002) (see also Figs. 5,7 and 10) and only in the northern portion of the Ejina Oasis.On the contrary, no gradient was observed in the wetter years (2000 and 2004).
The gradient in the northern portion of the Eijina Oases indicates a delay of about 20 days in the time of peak green-up over a distance of 15 km.Given the results in Fig. 11, this is most likely due to a combined influence of insufficient surface and ground water in the terminal section of the river and a gradient in temperature sum, with increasingly higher values occurring towards the river terminal zone.A map of the temperature sum, however, cannot be constructed due to the absence of air temperature observations at different locations in the Ejina Oasis.

Discussion
As shown in Sect.4, the annual vegetation growth is correlated with river streamflow upstream of the Ejina Oasis, both in the previous year and during the spring of the current year.The implication of such findings is that the shallow ground water is likely to be the major contributor of water resource to the vegetation growth in this area.
Even though the streamflow decreased since 1999 (Fig. 5), the water yield from the upper stream in 1998 was rather rich and local precipitation in the downstream area in 1999 was far above the long-term average.In particular, the positive anomalies in these two years occurred in the second half of the year when vegetation does not extract soil water, which most likely allowed ground water to be replenished sufficiently.Also, the temperature sum was the highest observed in the five years in 2000.As a consequence, vegetation conditions in 2000 were rather good, although both precipitation and the streamflow showed very large negative anomalies in 2000.It should be noted that the temperature sum is obviously determined by the current year evolution of weather conditions, so that the phenology of vegetation is determined by previous year hydrological conditions (most likely recharge of the shallow aquifer system) in Continuous rather low water yield from the upper stream in 1999 and 2000 (Fig. 5) and large negative anomaly for two years in local precipitation (2000 and 2001) (Fig. 4) resulted in 2001 in the worst vegetation growth condition during the period 2000-2005 as shown by the smallest annual mean NDVI and the delayed timing of peak green-up (Fig. 7).Water yield from the upper stream in 2002 was larger thanks to the new water allocation plan launched by the local government and as a consequence, a significant improvement was observed in the vegetation conditions in the following year 2003 (Fig. 8).Intensive water releases from the upper reach to the lower reach of the river in 2003 and moderate precipitation in 2003 resulted in fair good vegetation growth in 2004, although both the streamflow from upper reach and local precipitation showed a negative anomaly in 2004.
The combined effect of previous year water supply and of current year temperature sum is quite evident when considering the years 2000, 2001and 2002 (see Figs. 7, 10a and 11) (see Figs. 7, 10a and 11).Water supply was relatively limited during these years, with previous year streamflow being larger in 2000, which was characterized by warmer weather as indicated by the higher temperature sum.The period of vegetative growth was shorter, as indicated by the earlier peak green up, but vegetative development was better than in both 2001 and 2002, because of more abundant water.
Conversely, when comparing 2000 with 2003 and 2004, we see that the effect of previous year streamflow was dominant: weather was significantly cooler, i.e. the temperature sum significantly lower, leading to a much shorter duration of vegetative growth, but overall vegetative development was significantly higher than average due to the positive anomaly in streamflow (Fig. 10a).
According to hydrological studies by Su et al. (2007) and Wen et al. (2005), the ground water in the down stream of the Heihe River is recharged mainly from the rapid seepage flow of river water.Observations have shown that the time of propagation of fluctuations in river flow is determined by the distance to the river: up to 1 km from the river, where the shallow aquifer is recharged by current year's river flow.At distances between 1km to 10km from the river, the ground water recharge is a combination of the historical recharge (since 1950s) and recent years' recharge (Wen et al., 2005).These considerations, however, apply to the vicinity of each drainage channel.If we take the entire drainage system (see Fig. 1 in Akiyama et al., 2007) into account, we will have a combination of rapid propagation downstream along the drainage channels and slower outwards propagation from each drainage channel into the shallow aquifer.This mechanism is likely to explain the absence of time lag in peak green-up during the wetter years and the fact that a time lag does appear only in the northern, drier, portion of the Ejina Oasis and only in drier years (Fig. 13).
The depth of the shallow groundwater table fluctuates in response to fluctuations in streamflow, while dampening such fluctuations, which explains the absence of rapid fluctuations in the time series of the area-averaged NDVI.What needs to be evaluated in a future study is whether the amplitude of the fluctuations in groundwater table depth is as large as to have a significant impact on water uptake by roots and, therefore, on magnitude and timing of vegetation green-up.We can estimate the magnitude of the fluctuations in groundwater table depth on the basis of the observations of Akiyama et al. (2007)  groundwater table depth daily for 14 months at three sites 20 m (Site 1), 300 m (Site 2) and 10 km (Site 3) from the river.Yearly amplitude of shallow groundwater table depth was slightly less than 2 m at Site 1, less than 1m at Site 2 and nearly zero at Site 3, indicating fast recharge due to river water to at least 300 m from the river.Groundwater table depth was about 4 m at both Site 2 and Site 3. In Akiyama's study, estimated time of propagation of fluctuations in riverflow was 60 days for a distance of 300 m.The time lag estimated from the observations is about one month between Site 1 and 2 and several months between Site 1 and Site 3, although the very small seasonality at Site 3 makes this estimate less reliable.River water was observed occasionally at the river terminus during the non-irrigation months, while river water was present at the southern apex of the lower reaches during the non-irrigation periods.
Field observations in the Tarim Basin by Chen et al. (2004), Wang et al. (2007) showed non-stressed condition with water table at 3.82 m, mild stress with 4.74 m water table and moderate stress with 5.82 m.Similar studies in the HeiHe indicated a shallower optimal groundwater table depth, in the range 2 m to 4 m (Liu et al., 2007).
We may then conclude that wherever the groundwater table depth drops from 2m to 3 m say, a transition from wet, i.e. non limiting water supply, to dry, i.e. limiting water supply, conditions occur and we may expect to observe spatial variability in both magnitude and timing of vegetation green-up.Given the morphological and hydrological characteristics of the Ejina Oasis, this is more likely to occur in the northern portion of the Oasis.We may then tentatively explain the occurrence of the S-N gradient in the timing of vegetation green up by a decrease, say by 1 m, in groundwater table depth in the northern portion during drier years, so as to reduce significantly plant transpiration and photosynthesis.
The observed gradient in the 12-month phase value in our study (Fig. 13) implies that during drier years the green-up vegetation in the northern zone of the Ejina Oasis moves in S-N direction.This is consistent with the propagation of fluctuations in water supply: fast over longer distances through the drainage channels and slow over short distances through the shallow aquifer system.The literature evidence reviewed above suggest that "drier conditions", given the characteristics of the Ejina Oasis, means a decrease in groundwater table depth from 2m to 4m, possibly just 3m.Accordingly, our hypothesis, to be verified by a dedicate observations program, is that during years with lower streamflow groundwater table depth reaches to about 3 m in the northern portion but not in the southern portion, thus explaining the clear S-N gradient in the timing of vegetation green-up.

Conclusions
This study provided new insights on the hydrological processes determining vegetation phenology in the Ejina Oasis by looking at time series of ground based observations of precipitaion, air temperature and streamflow and of satellite observations of vegetation phenology.The primary driver of vegetation development is previous year streamflow, modulated by streamflow in the spring of the current year.Both overall development of vegetation, as measured by the yearly mean NDVI value, and the timing of peak green up are affected, but weather in the current year has a significant impact on both the duration and magnitude of vegetative development.This conclusion is based on the observed temperature sum and the time lag estimated by Fourier analysis of the NDVI time series.

4179
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | e. the Fast Fourier Transformation (FFT) (Menenti 4181 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | et al.

4183
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | If a time series of N samples of a variable Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | a ki is the original NDVI value in the k−th composite (i =1, 2, . . . . . .n) of the i −th year (k=1, 2, . . . . . .m ) in a pixel.HANTS algorithm uses all the available composite samples in the same year (along horizontal axis) to reconstruct time series of cloud-free NDVI data for this year.The TSS method focus on both temporal coordinates and uses the variance of the composite samples 4187 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | over different years to represent the similarity of the time series between the target year (with large gaps) and historical years: Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | ) 4189 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

4191
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | represents the annual mean value of NDVI.Analysis was focused on the vegetated land surface only by identifying pixels with land cover classes of forest, shrubs, sparsely Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | area-averaged annual mean NDVI values for the five years between 2000 and 2004 were analyzed against the annual streamflow of the pervious-2nd-year, the previous year and current year at the Zhengyixia hydrological station.A strong correlation 4193 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper |

4195
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

4197
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | combination with current year weather, most notably the temperature sum.

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper | in the Ejina Oasis.They measured 4199 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

4201
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | the reach of a large fraction of plant roots, thus leading to water stress and reduced growth.

Figure 2
Figure 2 Land cover map of the Ejina Oasis in 2004 (data provided by the Environmental and Ecological Science Data Center for West China, National Natural Science Foundation of China (http://westdc.westgis.ac.cn)).

Fig. 2 .
Fig. 2. Land cover map of the Ejina Oasis in 2004 (data provided by the Environmental and Ecological Science Data Center for West China, National Natural Science Foundation of China (http://westdc.westgis.ac.cn)).

Figure 3
Figure 3The temporal coordinate system for NDVI time series in a pixel.The horizonta for the composite coordinate in one year, while the vertical axis for year coordinate composite.ki a is the original NDVI value in the th k -composite ( k =1, 2, ……n) of th

Fig. 3 .
Fig.3.The temporal coordinate system for NDVI time series in a pixel.The horizontal axis stands for the composite coordinate in one year, while the vertical axis for year coordinate of the same composite.aki is the original NDVI value in the k-th composite (k=1, 2, ...n) of the i-th year (i=1, 2, ...m) at in a pixel.

Figure 4
Figure 4 (a) The annual precipitation in Dingxin station and Ejin Banner stations from 1980 to 2007; (b) The anomaly of the annual precipitations in Dingxin station and Ejina Banner station from 1980 to 2007.

Figure 4
Figure 4 (a) The annual precipitation in Dingxin station and Ejin Banner stations from 1980 to 2007; (b) The anomaly of the annual precipitations in Dingxin station and Ejina Banner station from 1980 to 2007.

Fig. 4 .
Fig. 4. (a) The annual precipitation in Dingxin station and Ejin Banner stations from 1980 to 2007; (b).The anomaly of the annual precipitations in Dingxin station and Ejina Banner station from 1980 to 2007.

Figure 5
Figure 5 (a) the annual stream flow and (b) the anomaly of the annual stream flow at Yingluoxia and Zhengyixia from 1980 to 2004.

Figure 5
Figure 5 (a) the annual stream flow and (b) the anomaly of the annual stream flow at Yingluoxia and Zhengyixia from 1980 to 2004.

Fig. 5 .
Fig. 5. (a) the annual stream flow and (b) the anomaly of the annual stream flow at Yingluoxia and Zhengyixia from 1980 to 2004.

Figure 6
Figure 6 (a) The monthly stream flow and (b) the anomaly of the monthly stream flow at Zhengyixia station from 1998 to 2004 (mean value was calculated between 1980 and 2004).

Figure 6
Figure 6 (a) The monthly stream flow and (b) the anomaly of the monthly stream flow at Zhengyixia station from 1998 to 2004 (mean value was calculated between 1980 and 2004).

Fig. 6 .Figure 7
Fig. 6.(a) The monthly stream flow and (b) the anomaly of the monthly stream flow at Zhengyixia station from 1998 to 2004 (mean value was calculated between 1980 and 2004).

Fig. 7 .Figure 8
Fig. 7. Area-averaged annual mean NDVI and phase value of 12-month Fourier harmonic component between 2000 and 2005 in the Ejina oasis of the Heihe River's down reaches.The phase value indicates the time of peak green-up of the vegetation.

Figure 9
Figure 9 (a) Area-averaged annual mean NDVI of year N versus total river flow of year (N-1); (b) Area-averaged annual mean NDVI of year n versus total river flow of March through May of year N.

Fig. 9 .
Fig. 9. (a) Area-averaged annual mean NDVI of year N versus total river flow of year (N-1); (b) Area-averaged annual mean NDVI of year n versus total river flow of March through May of year N.

Figure 10 Figure 10
Figure 10 (a) Area-averaged phase value of 12-month component of year N versus total river stream flow of year (N-1); (b) Area-averaged phase value of 12-month component of year N versus total river stream flow of April through May of year N. Streamflow was measured at the Zhengyixia gauge station.

Fig. 10 .Figure 11
Fig. 10.(a) Area-averaged phase value of 12-month component of year N versus total river stream flow of year (N-1); (b) Area-averaged phase value of 12-month component of year N versus total river stream flow of April through May of year N. Streamflow was measured at the Zhengyixia gauge station.

Fig. 11 .
Fig. 11.Area-averaged phase value of 12-month component of year N versus degree-days for May through July between 2000 and 2004.

Figure 12
Figure 12 The IHS color transformed image of Fourier harmonic components of NDVI from 2000 to 2008 in the Ejina Oasis in the down reaches of the Heihe river basin.'Intensity' = annual mean NDVI, 'Hue' = the phase of the 12-month harmonic component, 'Saturation' = the amplitude of the 12-month harmonic component.

Fig. 12 .Figure 13
Fig. 12.The IHS color transformed image of Fourier harmonic components of NDVI from 2000 to 2008 in the Ejina Oasis in the down reaches of the Heihe River Basin."Intensity" = annual mean NDVI, "Hue" = the phase of the 12-month harmonic component, "Saturation" = the amplitude of the 12-month harmonic component.

Fig. 13 .
Fig. 13.Column-average of phase values of 12-month component of year n versus row (i.e.S-N transect of oasis for each year N; gradient develops only on drier years; relate to fan-like drainage system and intermittently active drainage channels.
3425.98 E, and the bottom-right corner 41 • 50 6.85N,101• 1627.55E, which is shown in Fig.1.Ejina Oasis is the largest oasis in the lower reach of the Heihe River, characterized by a typical continental arid climate, with annual precipitation less than 100 mm and Inner Mongolia Provinces.The top-left corner of the study area is 42 • 16 46.64N, 101 • Zhu et al. (2009)estimate the critical range in shallow groundwater table depth.Zhu et al. (2009)studied the impact of fluctuations in groundwater table depth on water uptake by roots by means of numerical experiments performed with the model HYDRUS-1D and field experiments.In 2000 with average groundwater table depth at 2.64 m P.-euphratica, the dominant natural tree species in the area, obtained 53% of transpiration from groundwater during the central part of the growing season (DOY 160-290).The root vertical distribution was 10% in the top 60cm, 60% between 60 cm and 1 m, 25% between 1 m and 1.8 m and 5% between 1.8 m and 2 m.Numerical experiments indicated that a deeper water table depth from 2 m to 3 m would reduce transpiration by 74%.Deepening from 2.5 m to 3 m reduces transpiration from 1.5 mm d −1 to 0.5 mm d −1 during the central part of the growing season.Water flux from the groundwater table into the root zone was 300 mm with the groundwater table at 2 m, 150 mm with 2.5m and 60 mm with 3m.