Inter-annual and seasonal trends of vegetation condition in the Upper Blue Nile ( Abbay ) basin : dual scale time series analysis

Introduction Conclusions References


Introduction
Land degradation is a widespread environmental and development challenge (e.g.Dregne et al., 1991;UNEP, 2007).It is central to many international conventions and protocols related to environmental protection.Increasing demands for food, water and energy resulting from the growth in population and per capita consumption are driving unprecedented land use change (Godfray et al., 2010;Kearney, 2010).In turn, unsustainable land use is causing degradation of land resources.Thus, up-to-date quantitative information about land degradation is crucial to develop sustainable land use strategies and to support policy development for food and water security and environmental integrity.The status and trend of vegetation condition generally serve as a proxy for land degradation (Metternicht et al., 2010;Wessels et al., 2004Wessels et al., , 2007)).
Detecting and characterizing trends in vegetation condition over time using remotely sensed data has received considerable attention in recent years (Bai et al., 2008;de Jong et al., 2011;Eastman et al., 2013;Verbesselt et al., 2010a).Recent interest in vegetation trend analysis arises for three reasons.First, there is considerable interest in monitoring and assessing the state and trend of land degradation as well as in monitoring the performance of management programmes (Buenemann et al., 2011;Vogt et al., 2011).Second, it is only recently that substantial amounts of remotely sensed data and robust geospatial approaches suitable to such analysis have become available (Bai et al., 2008;Buenemann et al., 2011;Tucker et al., 2005).Third, it is a natural first step towards identifying drivers of changes in terrestrial ecosystems as vegetation variability and trends affect the exchange of water, energy, nutrients and carbon between the biosphere, the geosphere and the atmosphere (Baldocchi et al., 2001).
Changes in vegetation occur in three ways: (1) a seasonal or cyclic change that is driven by climate (e.g. annual temperature and rainfall) impacting plant phenology; (2) a gradual change over time that is consistent in direction (monotonic) such as change in land management or land degradation; and (3) an abrupt shift at a specific point in time (step trend) that may be caused by disturbances such as a sudden change in land use policies, deforestation, floods, droughts, and fires (Angert et al., 2005;de Jong et al., 2013a;Slayback et al., 2003;Tucker et al., 2001;Verbesselt et al., 2010b).Thus, consideration of all types of changes (i.e.seasonal, gradual and abrupt changes) is essential in order to assess the environmental impact of vegetation changes or to be able to attribute the changes in vegetation to drivers behind (de Jong et al., 2013a;Verbesselt et al., 2010b).There is substantial interest in monitoring the trends in seasonality of vegetation due to its sensitive response to climate change (Parmesan and Yohe, 2003;Sparks et al., 2009;Walther et al., 2002).
Long-term, remotely sensed normalized difference vegetation index (NDVI) data are suitable to detect and characterize trends in vegetation condition over time (de Jong et al., 2011;Eastman et al., 2013;Tucker et al., 2001;Verbesselt et al., 2010a).Although NDVI can be computed from different multispectral satellite data, NDVI from the Advanced Very High Resolution Radiometer (AVHRR) is the only global vegetation data set which spans a time period of three decades.Hence, it allows quantification of ecosystem changes as a result of ecosystem dynamics and varying climate conditions.In this study, AVHRR based Global Inventory, Monitoring, and Modeling Studies (GIMMS) NDVI data sets  were used for the purpose of long-term trend analysis (Tucker et al., 2005).In addition, NDVI time series (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011) from the Moderate Resolution Imaging Spectroradiometer (MODIS) 250 m (MOD13Q1) on board the Earth Observing System (EOS) Terra platform were used for medium-scale vegetation trend analysis.
Several vegetation trend studies from analysis of satellite observations have reported a positive trend in vegetation condition (i.e.greenness) in the Northern Hemisphere, including the Sahel (Fensholt et al., 2009;Karlsen et al., 2007;Slayback et al., 2003;Tucker et al., 2001).Other studies based on phenology and temperature observations from stations have also confirmed a greenness trend (Sparks et al., 2009).However, Zhao and Running (2010) reported a decreasing trend in greenness globally during the period 2000-2009 using Moderate-resolution Imaging Spectroradiometer (MODIS) NDVI data.Similarly, de Jong et al. (2011) found that many forested biomes experienced a decline in vegetation greenness.There is still no consistent result on trends of global vegetation activity (Bala et al., 2013).
The inconsistent results from various studies potentially emanated from differences in location of study areas (variations in altitude and latitude), trend detection techniques (ordinary least squares versus median trend), the length of the data series, and the NDVI data sources used (e.g.GIMMS, MODIS).Therefore, it is very important to characterize and understand inter-annual and intra-annual variability of vegetation activity using long-term data sets in study areas such as the Abbay Basin (1982Basin ( -2006)), where climate variability is high and topography is diverse.
It is essential to characterize and understand inter-annual and intra-annual variability of vegetation activity using longterm data sets such as GIMMS NDVI, especially if the increase or decrease in vegetation growth in the Abay Basin is mainly driven by climatic factors (De Beurs et al., 2009).However, if vegetation growth is mainly caused by human activities such as sustainable land management (SLM) programmes that involve localized implementation of best management practices (BMPs), then it is useful to use MODIS NDVI data as they provide repeated information at the spatial scale at which most human-driven land cover changes occur (De Beurs et al., 2009;Gallo et al., 2005;Townshend and Justice, 1988).Such analysis reveals areas of change that merit closer attention since it depicts significant shifts in local water, carbon and energy fluxes (Henebry, 2009;Vuichard et al., 2008).The results of this study are very im- portant from context of the Upper Blue Nile (UBN) Basin for two reasons.First, the UBN Basin is the major water source area of the Nile River, as it contributes around 70 % of the overall Nile flow.Thus, the UBN Basin is the most important river basin not only for Ethiopia but also for the downstream Nile Basin countries (i.e.Sudan and Egypt).Second, the UBN Basin accounts for a major share of the Ethiopia's irrigation and hydropower potential.
The general objective of this study was to characterize the degree of improvement or degradation of vegetation condition in the Abay Basin using both coarse-and medium-scale analysis.The specific objectives were (1) to identify interannual and seasonal trends in vegetation conditions at both coarse and medium scales using robust trend estimators and (2) to detect trend breaks.

Study area
The UBN/Abbay Basin covers a total area of 199 812 km 2 and is located in the centre and west of Ethiopia (Fig. 1).It lies approximately between 7 • 45 and 12 • 46 N, and 34 • 06 and 40 • 00 E. The altitude of the UBN Basin ranges from 475 m a.s.l. at the Sudanese border to 4257 m a.s.l. at the summit of Mt Guna.More than 83 % of the basin is located at an elevation above 1000 m a.s.l.The multi-year  mean annual rainfall varies between about 894 and 1909 mm a −1 , with a spatial mean of about 1396 mm a −1 based on gridded station-based rainfall data.The multi-year  mean annual areal potential evaporation ranges from 1065 to 1756 mm a −1 , with a spatial mean of about 1396 mm a −1 based on the Penman-Monteith method from the Climate Research Unit (CRU) 3.10 data set (Harris et al., 2014).The multi-year (1982The multi-year ( -2006) ) annual average temperature varies across space from 14 to 29 • C based on the CRU 3.10 data set.
Cropland is the major land use and land cover type, occupying more than 44 % of the total area of the basin (TEC- SULT, 2004).A variety of annual crops are grown under rainfed conditions, including wheat, barley, sorghum, teff, maize, finger millet, oil seeds, chick peas and beans.

GIMMS 15-day composite NDVI product
The Global Inventory, Monitoring, and Modeling Studies (GIMMS) NDVI data product (Tucker et al., 2005) was used to analyse the long-term  vegetation variability and trends in the UBN Basin.A release of the global coverage GIMMS data (GIMMSg) at 8 km resolution, covering from July 1981 through December 2006, was made freely available from the Global Land Cover Facility of the University of Maryland (http://www.landcover.org).The GIMMSg data set is composited at a 15-day time step by applying maximum value compositing (MVC) technique (Holben, 1986) in order to reduce cloud cover and water vapour effects.

MODIS NDVI 16-day L3 global 250 m (MOD13Q1)
The MODIS 16-day composite NDVI data sets with 250 m spatial resolution (MOD13Q1) tiles of h21v07 and h21v08 for the period 2001-2011 obtained from the NASA EOS data gateway were used in this study (Huete et al., 2002).The MODIS NDVI complements NOAA's AVHRR-derived NDVI products and provides continuity for time series historical applications.Apart from the vegetation index, the MOD13Q1 data contain pixel reliability information describing overall quality of NDVI values for each pixel.Pixel reliability information is a simple decimal number that ranks the product into five categories: −1, fill/no data; 0, good data; 1, marginal data; 2, snow/ice; and 3, cloudy.In the time series analysis, only "good" = 0 data are retained, and all other values are interpolated or replaced.

Land use and land cover (LULC) data
The third data set used in the analysis was LULC data produced by Woody Biomass Inventory and Strategic Planning Project (WBISPP) at a scale of 1 : 250 000 (TECSULT, 2004).The intention was to provide a basis for understanding the spatial distribution of seasonal trends.

Harmonic analysis of NDVI time series
Both AVHRR and MODIS NDVIs are composite products.Since maximum value composite (MVC) is not an atmospheric-correction method, some artefacts related to residual cloud cover and atmospheric haze remain in the MVC processed NDVI series (Nagol et al., 2009).Therefore, it was essential to perform harmonic analysis of both GIMMSg and MODIS NDVI data sets using the HANTS (Harmonic ANalysis of Time Series) algorithm (Roerink et al., 2000;Verhoef et al., 1996) for the purpose of (i) screening and removal of cloud contaminated observations and (ii) temporal interpolation of the remaining observations to reconstruct gapless images at a prescribed time.Harmonic (Fourier) analysis has been successfully applied to describing precipitation patterns (e.g.Landin and Bosart, 1989), meteorology (e.g.Legates and Willmott, 1990) and seasonal and interannual variations in land surface condition (Jakubauskas et al., 2001).The HANTS algorithm was developed based on the concept of discrete Fourier transform.HANTS considers only the most significant frequencies expected to be present in the time profiles, and applies a least-squares curvefitting procedure based on harmonic components (sine and cosine).Thus, harmonic (Fourier) analysis decomposes a time-dependent periodic phenomenon (e.g.vegetation development) into a series of sinusoidal functions, each defined by unique amplitude (strength) and phase (orientation with respect to time) (Roerink et al., 2000;Verhoef et al., 1996) and can be described as the sum of sine and cosine components as follows (Harris and Stöcker, 1998;Wilks, 2011): where f (t) is the series value at time t, T is the length of the series (e.g.T = 12 observations for a monthly time series), n is the number of harmonics to be used in the regression (in this study two harmonics were used), α 0 is the mean of the series, C n is the amplitude of harmonics given by C n = a 2 n + b 2 n , ω n is the frequency (2π n/T ), ϕ n is the phase angle given by ϕ n = tan −1 (b n /a n ) and e is the error term.Equation (1a) emphasizes the harmonic interpretation, while Eq. ( 1b) is a convenient transformation to a multiple linear harmonic regression model with coefficients a n = C n cos(ϕ n ) and b n = C n sin(ϕ n ) that can be easily estimated.
In this study the Interactive Data Language (IDL) implementation of the HANTS algorithm (De Wit and Su, 2005) was used.In HANTS algorithm a curve fitting is applied iteratively and the maximum iteration was set to 12 (Table 1).First, a least-squares curve is computed based on all data points, and next the observations are compared to the curve.Observations that are clearly below the curve are candidates for rejection, and the points that have the greatest negative deviation from the curve are therefore removed first.Second, a new curve is computed based on the remaining points and the process is repeated.Pronounced negative outliers are removed by assigning a weight of zero to them, and a new curve is computed.This iteration eventually leads to a smooth curve that approaches the upper envelope over the data points.Essential HANTS parameters were set as in Table 1 in order to obtain a reliable fitting curve.A detailed explanation of the parameters can be found in Roerink et al. (2000).Although the output comprises five Fourier com-  0, 1, 2 0, 24, 48  0, 1, 2 0, 23, 46  Invalid data rejection threshold: 0-1  0-1  0-1  0- ponents (i.e. three amplitudes and two phases), the zero frequency (mean annual NDVI) and the frequencies with time periods of 1 year (annual) were selected for the analysis.Amplitude 2 and phase 2 represent the semi-annual curve, but they are difficult to interpret (Eastman et al., 2013).This suite of harmonics (amplitude 0, amplitude 1 and phase 1) is sometimes called as the greenness parameters.

Long-term trend analysis
Before the time series analysis, the bimonthly HANTSfiltered GIMMS and MODIS data were aggregated to monthly composites by applying the arithmetic mean of each month.The bimonthly composites usually show strong serial autocorrelation, and aggregating into monthly observations reduces the effect of serial autocorrelation on the statistics generated from them.Seasonal variation must be removed in order to better discern the long-term trend in NDVI over time.Thus, a seasonally adjusted series is needed in order to better discern any trend that might have otherwise be masked by seasonal variation.Standardized anomalies (Z score) of the monthly HANTS-filtered NDVI data were calculated to remove the seasonality from the original time series (i.e.deseasoning) using the following equation: where x is the data value of the respective month, µ is the monthly long-term average value (i.e. the long-term average of all Januaries, Februaries, etc.) and σ is the monthly standard deviation.A Z-score value of 0 represents the long-term mean; a value of +1 represents 1 SD (standard deviation) of the long-term mean, and so on.The de-seasoned monthly NDVI time series was corrected for serial correlation (autocorrelation) using a trendpreserving Prewhitening (TPP) approach (Wang and Swail, 2001).The prewhitened series has the same trend as the original series, but with no serial correlation.Two types of interannual trend analysis and one trend significance test were applied to pre-whitened data on a pixel basis: linear trend (median trend or Theil-Sen), monotonic trend (Mann-Kendall) and Mann-Kendall significance test.
The median trend (Theil-Sen -TS) was computed for assessing the linear trend and rate of change per unit time (Hoaglin et al., 1983).The TS slope estimator was proposed by Thiel (1950) and modified by Sen (1968) and computes the median of the slopes between observation values at all pairwise time steps.The major advantage of the TS slope estimator is its breakdown bound.The breakdown bound for a robust statistic is the number of wild values that can occur within a series before it will be affected.For the median trend, the breakdown bound is approximately 29 %, meaning that the trends expressed in the image must have persisted for more than 29 % of the length of the series (in time steps) (Hoaglin et al., 1983).For example, in this study the GIMMS time series data contain 300 monthly images, so they are completely unaffected by wild and noisy values unless they persist for more than 87 months (i.e.0.29 × 300 = 7.25 years).The implication of this is that it also ignores the effects of short-term inter-annual climate teleconnections such as El Niño (typically a 12-month effect) and La Niña (typically a 12-24-month effect).

Seasonal trend analysis
Seasonal trend analysis was performed in order to identify trends in the essential character of the seasonal cycle while rejecting noise and short-term variability seasonal cycle (Eastman et al., 2013).This was performed based on a two-stage time series analysis.First, harmonic regression is applied to each year of images in the time series to extract an annual sequence of overall greenness and the amplitude and phase of annual and semi-annual cycles (Eastman et al., 2009).This suite of harmonics is termed as greenness parameters.Amplitude 0 represents the annual mean NDVI or overall greenness for each year.Amplitude 1 represents the peak of annual greenness.Phase 1 denotes the timing of annual peak greenness, represented by the position of the starting point of the representative sinewave of annual greenness.An increase or decrease in the phase angle means a shift in the timing to an earlier or later time of the year, respectively.The values of phase image potentially range from 0 to 359 • , such that each 30 • indicates a shift of approximately one calendar month (Eastman et al., 2009).Amplitude 2 and phase 2 represent the magnitude and position of a semi-annual curve, but their interpretation is difficult.They exist mainly as shape modifier of the annual curve (Eastman et al., 2013).Second, trends over years in the greenness parameters were analysed using a Theil-Sen median slope operator.This procedure is robust to short-term interannual variability up to a period of 29 % of the length of the series (Eastman et al., 2013;Gilbert, 1987;Hoaglin et al., 1983).Furthermore, the samples are annual measures of shape parameters.Thus, short inter-annual disturbances (i.e.those 8 years or less in this case) such as individual El Niño-Southern Oscillation (ENSO) events, have little effect on the long-term trends represented by the median trends.The significance of the trends in the five shape parameters was tested using the contextual Mann-Kendall (CMK) statistics test (Neeti and Eastman, 2011), which is a modified form of the Mann-Kendall test of trend significance (Kendall, 1948;Mann, 1945).The CMK trend test reduces the detection of spurious trends and an amplification of confidence when consistent trends are present through adding contextual information.CMK is a modified version of the Mann-Kendall test which is based on a principle that a pixel would not be expected to exhibit a radically different trend from neighbouring pixels.
In order to make the interpretation easier, the seasonal trends were categorized into different classes based on major CMK trend significance.The CMK trend significance for each of the five shape parameters was classified into three major categories: significantly decreasing at the p < 0.05 level, not significant, and significantly increasing at the p < 0.05 level.Only the major seasonal trend classes in terms of area were selected for detailed analysis.

Change-point analysis
A change-point analysis was conducted, because it indicates that the trends change between positive and negative within the analysis period.Regions of interest (ROIs) were selected based on most prevalent classes of significant seasonal trends for GIMMS NDVI and MODIS NDVI.The Breaks For Additive Seasonal and Trend (BFAST) method (Verbesselt et al., 2010a, b) was used to detect and characterize abrupt changes within the trend and seasonal components of the ROI time series.BFAST is particularly important if long-term trends are composed of consecutive segments with gradual change, separated from each other by a relatively short period of abrupt change (de Jong et al., 2012;Verbesselt et al., 2010b).In BFAST, the ordinary least-squares (OLS) residuals-based moving sum (MOSUM) test is used to test for whether one or more breakpoints are occurring.Thus, the use of BFAST trend-break analysis assisted identification of several change periods which might otherwise been overlooked through a commonly assumed fixed change trajectory analysis.BFAST is an iterative algorithm that integrates the decomposition of time series into seasonal, trend, and remainder component with methods for detecting changes.An additive decomposi-tion model is used to iteratively fit a piecewise linear trend and a seasonal model.The general model is of the form Y t = T t + S t +e t , (t ∈ m), where Y t is the observed NDVI data at time t (in the time series m), T t is the trend component, S t is the seasonal component, and e t is the remainder component.A trend component corresponds to a longterm process that operates over the time spanned by the series.A remainder component corresponds to a local process which causes variability between cycles.A seasonal component corresponds to a cyclic process that operates within each cycle.

GIMMS NDVI
The Upper Blue Nile (Abay) Basin is characterized by a multi-year mean NDVI of 0.49 for the period 1982-2006 using the GIMMS data set.On the one hand, the south and southwest regions of the basin represent the highest mean NDVI and the lowest temporal variation.This area is characterized by shrubland and woodland vegetation.On the other hand, the eastern, northeastern and northwestern parts represent the lowest mean NDVI and the highest temporal variation as a consequence of the relatively low annual precipitation and the more extended dry periods.Moreover, this area is characterized by intensive and sedentary cultivation.
During the period 1982-2006, a significant increasing trend in vegetation condition is observed in large part of the basin, particularly in the western and central parts of the basin (Fig. 2).About 77 % of the basin showed a positive trend in monthly NDVI with a mean rate of 0.018 NDVI yr −1 (i.e.0.0015 × 12 NDVI units), and the remaining 23 % of the basin exhibited a decreasing trend with a mean rate of 0.0084 NDVI yr −1 (i.e.0.0007 × 12) (Table 2).About 41 % of the basin area showed a significant (p < 0.05) increase in monthly NDVI, with a mean rate of 0.03 NDVI yr −1 (0.0023 × 12), and about 4.6 % of the basin showed a significant decrease in monthly NDVI, with a mean rate of 0.02 NDVI yr −1 (0.0017 × 12).

MODIS NDVI
Areas of higher multi-year mean NDVI depict lower vegetation cover variability (Fig. 3a and b).Areas located in the eastern part of the basin are those with higher variability and show an upward (positive) trend of vegetation growth in the period 2001-2011 (Fig. 3b and c).About 36 % of the UBN/Abay Basin shows a significant browning trend (p < 0.05) over the period 2001-2011 at an average rate of 0.0768 NDVI yr −1 using MODIS 250 m spatial resolution (Tables 3 and 4).By contrast, only about 1.19 % of the basin shows a significant greening trend (p < 0.05) over the period 2001-2011 at an average rate of 0.066 NDVI yr −1 (Ta-  bles 3 and 4).The majority of the browning trend comes from Dabus (14.27 % at a rate of 0.0744), Belles (13.5 % at a rate of 0.0816), Wonbera (11.81 % at a rate of 0.0792) and Didessa (11.48 % at a rate of 0.0792).The vegetation activity in the Belles sub-basin is decreasing at a faster rate (0.0816 NDVI yr −1 ) than the others (Tables 3 and 4).The maximum significant browning slope observed in the Didessa sub-basin is the largest of all maximum slopes in the basin.The size of greening and the browning area within sub-basins are proportional in sub-basins such as Beshillo, Jemma, Muger and Welaka.However, more than half of the areas of the Belles (66.91 %), Rahad (60.41 %), and Wonbera (64.12 %) are represented by a significant browning trend.The majority of the greening trend comes from the Jemma (22.03 % at a rate of 0.066), Beshillo (18.23 % at a rate of 0.066), North Gojam (22.78 % at a rate of 0.066), and Muger (12.15 % at a rate of 0.0648) sub-basins.The slowest and fastest rates of significant greening are observed in the Dinder sub-basin (0.0552 NDVI yr −1 ) and Dabus sub-basin (0.0768 NDVI yr −1 ), respectively.

GIMMS NDVI
More than half (59.52 %) of land areas exhibit a significant trend in seasonality over the 25-year time period measured , producing a total of 53 significant trending classes of seasonal curves which include different combinations of significant changes in amplitude 0, amplitude 1, amplitude 2, phase 1 and phase 2. Thus, a trend change in seasonality is not a rare occurrence in the Abay Basin.Out of the 53 classes, only the prominent 5 classes of seasonal trends are described below, which account for about 67 % of all significantly trending pixels and represent 39.95 % of all land pixels of the basin.Figure 4 depicts the spatial distribution of the five most prevalent classes of significant seasonal trends in GIMMS NDVI.
www.earth-syst-dynam.net/6/617/2015/Earth Syst.Dynam., 6, 617-636, 2015  Class 1 has an area of 19 712 km 2 , which accounts for 16.85 % of all significantly trending areas and 10.03 % of the total land area (Table 5).This class represents areas with a significantly positive trend in amplitude 0 (mean annual NDVI) and no significant change in any of the other four parameters (amplitudes, phase 1 and phase 2), producing curves that are raised uniformly throughout the year.The frequency distribution of seasonal trend classes within zones of land cover types shows that the majority of class 1-trending pixels occur in areas of woodland (33 %) and shrubland (30 %).Thus, it appears that woodland and shrubland areas are experiencing a consistent increase in productivity throughout the year.
Class 2 is characterized by both a significant increase in amplitude 0 (mean annual NDVI) and a significant decrease in phase 1 (orientation with time of the annual cycle).A de-Table 3. The percentage of a sub-basin's significant trend out of the total area of significant trend in the basin and the percentage of area significantly increasing or decreasing within sub-basins are shown.The result is based on MODIS-based analysis.The locations of the sub-basins can be found in Fig. 8  growing season and no effect on the start of the growing season, leading to lengthening of the growing season.Class 3 represents areas that exhibit a significant increase in amplitude 1 (the annual cycle) but no other significant trend parameters.Class 3 covers 3.48 % (6848 km 2 ) of all the basin area (Table 5).This class of seasonal trend exhibits a characteristic of increase in NDVI during the growing season (August-November) balanced by declines during the dry season (February-May) in order to maintain the same mean NDVI (Fig. 5c).The majority of class 3 trending pixels occur in cropland areas (68 % or 4672 km 2 ) of the basin.Therefore, NDVI is enhanced during the green season (cropping season) and inhibited during the brown season or harvest.
Class 4 represents areas that exhibit a significant decrease in amplitude 1 (the annual cycle) but no other significant trend parameters.Class 4 covers 9.51 % (18 688 km 2 ) of all the basin area (Table 5).This class of seasonal trend exhibits a characteristic decrease in NDVI during the growing season (more pronounced between August and November) increase in NDVI during dry season (Fig. 5d).Class 4 seasonal trend pixels are mostly associated with cropland (59 %) and grassland (19 %) land cover types.This type of seasonal trend is more common in the North Gojam and Lake Tana sub-basins.
Class 5 describes areas that show a significant decrease in phase 1 only.Class 5 seasonal curves resemble class 2 seasonal curves in that both have similar pattern of a decrease in yearly phase.It covers 7.5 % (14 720 km 2 ) of all the basin area and represents 12.6 % of all significant seasonal trend pixels.This class of seasonal trend exhibits a characteristic shift at the end of growing season.Class 5 pixels are mostly associated with woodland (48 %) and cropland (30 %).This type of seasonal trend is more common in the southern and southwestern parts of the basin.

MODIS NDVI
About 207 classes of significant intra-annual changes were obtained from the MODIS-based seasonal trend analysis during the period 2001-2011.These classes include different combinations of significant changes in greenness parameters.Approximately 42 % of the Abay Basin shows significant seasonal trends, out of which 96 % is comprised of the dominant eight classes (Table 6 and Fig. 6).Almost all vegetated land cover types commonly and dominantly exhibit the seasonal trend characteristics of class 1, class 2 and class 3 (Table 6).In addition, a class 4-type seasonal trend is observed mainly over 14 % of the woodland and 15 % of the Earth Syst.Dynam., 6, 617-636, 2015 www.earth-syst-dynam.net/6/617/2015/all significant pixels (p < 0.05) (Table 6).This class represents areas with a significant decreasing trend in peak of annual greenness (amplitude 1), with no significant changes in overall greenness and time peak greenness during the period 2001-2011.Figure 7a shows a decrease in greenness in the wet season (June, July, August, September) but an increase in greenness in the dry season (January, February, March, April).Class 2 represents significant decreases in overall greenness and annual peak greenness, with no significant changes in other greenness parameters (Table 6).Figure 7b shows the seasonal curves of NDVI (2-year average) for class 2. The maximum difference between the two curves occurs in August, whereas the minimum difference occurs in February and November.NDVI values for both curves are above 0.6.Class 3 has an area of 11 408 km 2 and describes 16.21 % of all significant seasonal trend and 6.8 % of the basin total area (Table 6).Class 3 represents areas with a significantly negative trend in overall greenness (amplitude 0) and no significant change in any of the other greenness parameters, producing a curve that decreases uniformly throughout the year.The NDVI values for all months remain above 0.7 throughout the year.About 34 % of the natural forest corresponds to class 3 type of seasonal trend.Thus, the change in overall greenness may indicate a decrease in vegetation density without complete removal forest.Class 4 describes about 8 % of all significant trends and covers an area of 6525 km 2 (Table 6).Class 4 includes areas that exhibit a significant increase in phase 1 (the time of annual peak greenness) but no other significant seasonal trend parameters.Figure 7d   creasing or decreasing but can change over time.For example, ROI #2 (Fig. 8b) and ROI #4 (Fig. 8d) represent gradual browning and greening trends, respectively, over the entire study period, while ROI #1 and ROI #3 describe abrupt browning and greening, respectively.ROI #1 is characterized by trend breaks at December 1989 and August 1997 and ROI #3 shows a trend break at April 1998.Thus, the trend-break analysis assisted identification of several change periods which might be overlooked through a commonly assumed fixed change trajectory analysis.

Trend-break analysis in the HANTS seasonality parameters
Trend-break analysis was performed on the harmonic shape parameters to identify the time when abrupt change occurred.No significant break points were detected based on the OLSbased MOSUM test for structural change using BFAST.However, through visual interpretation it is possible to identify trend-break points.For example, class 2 and class 5 showed a sharp increase in phase 1 between 1988 and 1998 and a sharp decrease afterwards (Fig. 9).

MODIS NDVI
The trend-break analysis did not indicate significant trend breaks using the BFAST approach in all of the sampled areas (Fig. 10), suggesting that vegetation activity in the Abay Basin showed a monotonic increasing or decreasing behaviour.Slowly acting processes such as land management practices or land degradation may cause such monotonic changes in the time series.ROI #1 and ROI #2 are similar in representing areas of significant monotonic decreasing trend in vegetation condition but different in the reasons behind the trend.ROI #1 is located in the Didessa sub-basin, and the main reason for the observed decreasing trend is the decrease in woody vegetation density as a result of human activity.ROI #2 is located in the Beless sub-basin and forest fire is the main cause of the observed decrease in vegetation condition.Both ROI #3 and ROI #4 represent areas of significant increasing trend in veg-

Trend-break analysis in the HANTS seasonality parameters
Although there is a significant monotonic increasing or decreasing trend, no significant break points were identified in any of the greenness parameters based on the OLS-based MOSUM test.The trends in overall greenness and annual peak greenness appear to be continuously increasing or decreasing without trend breaks.However, phase 1 in class 5, class 6 and class 8 shows a non-significant trend break (Fig. 11).The time of annual peak greenness has increased since the year 2007.

GIMMS-based inter-annual and seasonal trends
Vegetation trend analysis using the monthly GIMMS NDVI revealed a significant increase in the vegetation condition over 41.15 % of the Abay Basin for the period 1982-2006.The positive trend in the eastern part of the basin could be related to the rehabilitation effort on degraded areas by the government as a response to the devastating drought in 1984/1985.The northwestern and western parts of the basin that exhibited positive NDVI change are mainly covered with lowland vegetation such as woodland and shrubland.These areas are used for traditional purposes such as shifting cultivation and livestock rearing by the small indigenous population.One key characteristic of these areas is the periodic burning of vegetation, as is evident from fire scars visible on Landsat images.The practice of burning stimulates the growth of grass in the burnt area, and it can provide forage to be grazed for the next season.Thus, the positive NDVI trend of these areas could be related to the decrease in the trend www.earth-syst-dynam.net/6/617/2015/Earth Syst.Dynam., 6, 617-636, 2015 of burning of woody vegetation.Several studies have also reported widespread greening and increased vegetation growth and productivity for various areas of the world, particularly in the Northern Hemisphere (Nemani et al., 2003;Slayback et al., 2003;Zhou et al., 2001).However, detailed investigation of these areas is crucial in order to be able to attribute the trend to particular cause.The basin is characterized not only by significant interannual variation but also by a significant trend in seasonality.More than half (59.52 %) of land areas exhibited a significant trend in seasonality over the 25-year time period measured .Thus, a trend in seasonality is not a rare occurrence in the Abay Basin.In their global study, Eastman et al. (2013) also concluded that a significant trend in seasonality is not a rare occurrence.Areas with a significantly positive trend in mean annual NDVI (amplitude 0) were found to be the dominant seasonal trend class, which accounts for 16.85 % of all significantly trending areas and 10.03 % of the total land area.The majority of this seasonal trend class occurs in areas of woodland and shrubland.A significant increase in mean annual NDVI without significant changes in any of the other harmonic shape parameters, as amplitude 0 correlates with annual integrated NDVI rescaled to a mean value.Integrated NDVI correlates strongly with annual gross primary productivity (Reed et al., 1994).Thus an increasing trend in amplitude 0 would imply increased productivity uniformly throughout the year.It appears that woodland and shrubland areas of the Abay Basin are experiencing a consistent increase in productivity throughout the year.This result is consistent with previous findings over the African Sahel region (Eklundh and Olsson, 2003;Anyamba and Tucker, 2005;Olsson et al., 2005;Heumann et al., 2007) The most likely driving factors of monotonic greening or browning could be changes in growth-limiting climatologies (Nemani et al., 2003), broad-scale land management practices, and persistent land use change driven by settlements (e.g.cropland expansion and urbanization) (Ramankutty et al., 2007).The strongest indication for greening in agricultural expansion areas such as the Fincha and Belles subbasins is most likely attributable to improved agricultural techniques.Abrupt browning could be caused by logging followed by regrowth or the already existing vegetation decline in some places might have amplified by a period of persistently poor rains.However, abrupt greening could be caused by a wet period induced by El Niño/La Niña-Southern Oscillation warming events such as ENSO 1986/87, 1994/95, and 1997/98. The year 1997/1998 was identified as a trendbreak point using the Pettitt test.Thus, greening might represent recovery from drought or other disturbances such as forest fire (Anyamba and Tucker, 2005;Heumann et al., 2007;Olsson et al., 2005).Greening as a result of land management such as forest plantation cannot be traced with GIMMS data, as most of the plantations are situated in the area at the local scale and individual farm level.Finally, although the GIMMS data sets have been thoroughly corrected (Tucker et al., 2005), AVHRR satellite platform changes could potentially cause trend breaks within the time series.The detected trend-break time was compared with the time of AVHRR platform change (1985, 1988, 1994, 1995, 2000 and 2004) (Tucker et al., 2005).With regard to attribution of vegetation trend to the drivers, the coarser-scale analysis only provides indications for more focused analysis, because the driving factors might act on a local scale much smaller than the resolution of GIMMS NDVI (i.e. 8 km).Thus, based on the information obtained from this courser-scale analysis, areaspecific attribution of the increase in greenness to its causes is required for the Abay Basin for the future.

MODIS-based inter-annual and seasonal trends
The vegetation trend analysis with a spatial resolution of 250 m MODIS (2001-2011) depicted a different perspective from that of 8 km spatial resolution of GIMMS data set .About 36 % of the UBN/Abay Basin showed a significant browning trend, and only 1.19 % of the basin showed a significant greening trend over the period 2001-2011.The majority of the browning trend comes from Dabus, Belles, Wonbera and Didessa.The MODIS-based vegetation seasonal trend analysis during the period 2001-2011 revealed four conspicuous vegetation changes, accounting for 86 % of all significant pixels: (1) a change in the greenness pattern from a strong annual cycle to strong semi-annual vegetation cycle, indicating significant changes from shrubland/woodland vegetation to double cropping (irrigation); (2) a decrease in both overall greenness and annual peak greenness; (3) degradation of vegetation without complete removal (e.g.ROI #3); and (4) a change from semi-annual vegetation cycle to strong annual vegetation cycle.These areas could have been masked if only the inter-annual vegetation trend analysis were performed.Thus, carrying out intraannual trend analysis at medium resolution (250 m) would be helpful since most of the deforestations for the purpose of small-scale irrigation are operating at the local scale.However, a mere identification of areas with significant intraannual vegetation change is not enough for monitoring of vegetation condition.Attribution of the observed categories of intra-annual changes to the driving factors based on detailed field observation is also important.
The BFAST approach of the trend-break analysis indicated a monotonic increasing or decreasing vegetation condition without any significant trend breaks.The greening trend observed from the GIMMS data set might have been reversed into a browning trend, suggesting a localized decline in vegetation activity during the period 2004-2006.This result is consistent with the findings of de Jong et al. (2013b), who explained a trend reversal into browning trend not confined to single geographical region but extended across all continents.In the 1980s and 1990s the Northern Hemisphere was found to become greener (Nemani et al., 2003;Zhou et al., 2001).The impact of sensor degradation on trends in MODIS NDVI (Wang et al., 2012) could be one cause of the recent browning trend in the basin.Other studies attributed anthropogenic factors to the observed recent browning trend in the basin.In the northwest part of the basin, especially in the Rahad, and some parts of Dinder, sub-basins could be affected by intra-regional resettlement (i.e.people from the same area with kin relations in the same locality) and its associated deforestation for cropland expansion (Dixon and Wood, 2007;Lemenih et al., 2014).The Amhara regional government of Ethiopia initiated a "voluntary" resettlement scheme in 2003 for the most chronically food-insecure people from all zones of the region to potentially more productive, fertile and less populated parts, such as the Metema, Quara, Tach Armacheho and Tegede weredas of North Gondar administrative zone (NCFSE, 2003).In the Belles sub-basin, land preparation for small-to medium-scale irrigation schemes could initially cause a decline in vegetation trend, but after a while the browning trend might be reversed to a greening trend.Currently, private and government agricultural investments are underway in this sub-basin.Thus, the decreasing trend of vegetation greenness in the Belles sub-basin can be largely explained by the de-vegetation during land preparation for sugarcane plantation, particularly in the Alefa and Jawi weredas of Amhara regional state and the Pawe and Dandur weredas of Benshangul-Gumuz regional state of Ethiopia.Such large croplands could contribute to the brown- Intra-annual trend analysis identified changes in vegetation condition that could have been masked if only the inter-annual vegetation trend analysis were performed.Intraannual trend in vegetation condition was found to be prevalent in the landscape of Upper Blue Nile/Abay Basin based on both the GIMMS-based and MODIS-based analysis.Only five types of seasonal trends were dominant, out of which the largest proportion describes areas that are experiencing a uniform increase in NDVI throughout the year.It appears that woodland and shrubland areas of the basin experienced an increase in vegetation productivity resulting from a longer growing season.
The trend-break analysis of GIMMS NDVI (1982NDVI ( -2006) ) showed that trends in vegetation conditions of the UBN Basin were found to be not only monotonic but also abrupt.However, the trend-break analyses performed on MODIS NDVI (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011) showed that vegetation condition in the basin is a monotonic increasing or decreasing trend.Slowly acting processes such as land management practices or land degradation may cause such monotonic changes in the time series.

Figure 1 .
Figure 1.Location map of the study area, Upper Blue Nile/Abay Basin, Ethiopia.
Eucalyptus plantation mainly for wood and fuel wood is common in areas near villages.Extensive areas of bamboo occur in the lower areas of the western part of the basin.Dense woodland is mainly found on the lower western slopes in Wellega, North Gondar and Benishangul-Gumuz; open woodland is found notably in Wellega, Illubabor and Jimma; and open shrubland is found in highland areas.Woodlands and shrublands are often associated with a grass under-storey or grassed areas between areas of low woody vegetation.Two main types of grassland occur in the basin: (a) lowland tall grasslands occur in low-rainfall areas interspersed with a few 620 E. Teferi et al.: Inter-annual and seasonal trends of vegetation condition in the Upper Blue Nile (Abay) Basin trees and shrubs, and (b) highland temperate grasslands occur mainly above 2000 m a.s.l.

Figure 2 .
Figure 2. The inter-annual trend in vegetation condition based on GIMMS NDVI in the Abay (UBN) Basin for the period 1982 to 2006.

Figure 3 .
Figure 3. Vegetation cover variability and trends for the Upper Blue Nile (Abay) Basin based on MODIS 250m NDVI (MOD13Q1) for the period 2001-2011: (a) multi-year mean monthly NDVI, (b) coefficient of variation of monthly NDVI, (c) Theil-Sen slope ( NDVI/month), and (d) Mann-Kendall significance test of the trend test.

Figure 4 .
Figure 4.The five most prevalent classes of significant seasonal trends in GIMMS NDVI: class 1 shows significant increases in amplitude 0 (mean NDVI), class 2 shows significant increases in amplitude 0 and decrease in phase 1 (the timing of annual peak greenness) together, class 3 shows significant increases in amplitude 1 (increase in the difference between minimum and maximum NDVI without affecting the mean), class 4 shows significant decreases in amplitude 1, and class 5 shows significant decreases in phase 1.

Figure 5 .
Figure 5. Seasonal curves representative of each of the major trend classes from the 1982-2006 series: the green curve represents the characteristic seasonal curve at the beginning of the series (mean value of 1982 and 1983), while the red curve represents the characteristic seasonal curve at the end of the series (mean value of 2005 and 2006).

Figure 6 .
Figure 6.The most prevalent classes of significant seasonal trends of monthly 250 m MODIS NDVI (MOD13Q1) for the Upper Blue Nile/Abay Basin from 2001 to 2011.

Figure 7 .
Figure 7. Fitted seasonal curves for 2001-2002 in green and 2010-2011 in red, derived from trends over the complete series.The intraannual variations of generalized monthly NDVI over the first and last 2-year periods of the time series (2001-2002 and 2010-2011) are shown for the eight ROIs selected based on the amplitude composite images (Fig. 8).

Figure 8 .
Figure 8. Examples of decomposition and trend-break analysis of monthly GIMMS NDVI time series for sampled regions of interest (ROIs) of increasing trend (ROI #1 and ROI #2) and decreasing trend (ROI #3 and ROI #4) generated by the BFAST approach.The top section in every plot shows the NDVI data, whereas the other three sections depict the individual components after decomposition.The seasonal (S t ) and remainder (e t ) components have zero mean, while the trend component (T t ) shows the trend in NDVI.The slope coefficients (β) and the significance levels (P ) at α value of 0.05 for each segment are given.
depicts the seasonal curves of class 4 during the initial (average of2001 and 2002)  and final period (average of 2010 and 2011) of the analysis.This class of seasonal trend exhibits a shift in the timing of annual peak greenness to an earlier time of the year.Class 5 represents areas with a significant decrease in both amplitude 0 and amplitude 1 but significant increase in phase 1.This type of seasonal trend class is commonly found in the Wonbera sub-basins.The NDVI values remain above 0.6 in all months throughout the year.According to Fig.7e, the maximum difference between curves occurs between June and November, suggesting that the density of land cover type causing this difference might have been reduced.Class 6 represents areas with a significantly negative trend in phase 1 and no significant change in any of the other greenness parameters.This class is mostly found in the North Gojam sub-basin and is characterized by lower NDVI values throughout the year as compared to the NDVI values of other classes.Class 7 exhibits a significant decrease in overall greenness and significant increase in the time peak annual greenness.The NDVI values for all months remain above 0.7 throughout the year, and this makes the seasonal curve of class 7 similar to that of class 3.Both class 3 and class 7 are mostly located in the Didessa sub-basin.Class 8 represents areas with significant decrease in annual peak greenness and significant increase in time of peak annual greenness.These changes are noticeable in Fig.7h.The annual peak greenness observed in the initial state(2001/2002)  decreased by approximately 0.2 NDVI units between July and October as compared to the final state (analysis in the inter-annual data seriesThe trend-break analysis of ROIs shown in Fig.8illustrates different trend behaviour within a longer time series of GIMMS data.Trends are not always continuously inwww.earth-syst-dynam.net/6/617/2015/Earth Syst.Dynam., 6, 617-636, 2015

Figure 9 .
Figure 9. GIMMS-based trends and trend changes on harmonic series shape parameters for (a) class 1 (significant increases in amplitude 0 -mean NDVI); (b) class 2 (significant increases in amplitude 0 and decrease in phase 1 (the timing of annual peak greenness) together); (c) class 3 (significant increases in amplitude 1); (d) class 4 (significant decreases in amplitude 1); and (e) class 5 (significant decreases in phase 1).No breakpoints are detected in all classes based on the ordinary least-squares (OLS) residuals-based moving sum (MOSUM) test in BFAST.

Figure 10 .
Figure 10.Examples of decomposition and trend-break analysis of monthly MODIS NDVI time series (2001-2011) for certain regions of interest (ROIs) for significant decreasing trends, (a) ROI #1 and (b) ROI #2, and for significant increasing trends, (c) ROI #3 and (d) ROI #4, generated by the BFAST approach.The locations of ROIs are shown in Fig. 3d.The top section in every plot shows the NDVI data, whereas the other three sections depict the individual components after decomposition.The seasonal (S t ) and remainder (e t ) components have zero mean, while the trend component (T t ) shows the trend in NDVI.The slope coefficients (β) and the significance levels (P ) at α value of 0.05 for each segment are given.

Figure 11 .
Figure 11.MODIS-based trends and trend changes on harmonic series shape parameters for most prevalent classes of significant seasonal trends.The sampled locations are displayed in Fig. 8.
www.earth-syst-dynam.net/6/617/2015/Earth Syst.Dynam., 6, 617-636, 2015 E. Teferi et al.: Inter-annual and seasonal trends of vegetation condition in the Upper Blue Nile (Abay) Basin ing in the initial stage and to the greening once the cropland is established.5 Conclusions This study concludes that integrated analysis of inter-annual and intra-annual NDVI trend based on GIMMS and MODIS data is crucial for robust identification of changes in vegetation condition.The GIMMS-based inter-annual trend analysis revealed that about 41 % of the UBN Basin depicts a significant increase in monthly NDVI during the period 1986-2006.However, the MODIS-based inter-annual trend analysis revealed that about 36 % of the UBN Basin shows a significant decreasing trend in NDVI over the period 2001-2011.The greening trend observed from the GIMMS data in the 1990s might have been reversed into a decreasing trend in the 2000s.Other studies in the area have attributed anthropogenic factors (e.g.resettlement and agricultural expansion) to the observed greening-to-browning reversal of vegetation conditions in the UBN Basin.The decreasing trend in vegetation condition observed from MODIS-based analysis could indicate a localized decline in vegetation activity in recent years.

Table 1 .
Parameters used in HANTS analysis.

Table 2 .
Result of linear trend analysis based on Theil-Sen (TS) (median) slope for the monthly GIMMS NDVI data.
Absolute change refers to the median slope; relative change = (median slope/multi-year mean) × 100.

Table 4 .
. Summary statistics of Theil-Sen slope (median trend) for the significantly increasing (greening) and decreasing (browning) trends per sub-basin.Locations of the sub-basins can be found in Fig.8.The minimum, maximum and mean of significantly browning slopes have negative signs.The result is based on MODIS-based analysis.creasingphaseanglemeansashift to a later time of the year.This class has an area of 18 560 km 2 and describes about 10 % of the basin land area (Table5).Figure5bshows that the green-down period is occurring about 15 days later in 2006 than 1982.About 56 % of class 1-trending pixels occur in areas of woodland.This means that woodland areas of the Upper Blue Nile/Abay Basin are experiencing an increase in mean annual NDVI with a significant shift at the end of the www.earth-syst-dynam.net/6/617/2015/Earth Syst.Dynam., 6, 617-636, 2015

Table 5 .
Percent of seasonal trend classes by land cover classes for both GIMMS-based and MODIS-based analysis.

Table 6 .
Extent of dominant classes of trends in seasonality from MODIS-based analysis.