Global changes in dryland vegetation dynamics (1988–2008) assessed by satellite remote sensing: comparing a new passive microwave vegetation density record with reflective greenness data

Drylands, covering nearly 30 % of the global land surface, are characterized by high climate variability and sensitivity to land management. Here, two satellite-observed vegetation products were used to study the long-term (1988– 2008) vegetation changes of global drylands: the widely used reflective-based Normalized Difference Vegetation Index (NDVI) and the recently developed passive-microwavebased Vegetation Optical Depth (VOD). The NDVI is sensitive to the chlorophyll concentrations in the canopy and the canopy cover fraction, while the VOD is sensitive to vegetation water content of both leafy and woody components. Therefore it can be expected that using both products helps to better characterize vegetation dynamics, particularly over regions with mixed herbaceous and woody vegetation. Linear regression analysis was performed between antecedent precipitation and observed NDVI and VOD independently to distinguish the contribution of climatic and non-climatic drivers in vegetation variations. Where possible, the contributions of fire, grazing, agriculture and CO 2 level to vegetation trends were assessed. The results suggest that NDVI is more sensitive to fluctuations in herbaceous vegetation, which primarily uses shallow soil water, whereas VOD is more sensitive to woody vegetation, which additionally can exploit deeper water stores. Globally, evidence is found for woody encroachment over drylands. In the arid drylands, woody encroachment appears to be at the expense of herbaceous vegetation and a global driver is interpreted. Trends in semi-arid drylands vary widely between regions, suggesting that local rather than global drivers caused most of the vegetation response. In savannas, besides precipitation, fire regime plays an important role in shaping trends. Our results demonstrate that NDVI and VOD provide complementary information and allow new insights into dryland vegetation dynamics.


Introduction
Drylands cover nearly 30 % of the global land surface, they are characterized by high climate variability and are sensitive to land use practice (Tietjen et al., 2010). Over the last decades, many dryland ecosystems have faced increased pressure from human demands and climate change (Asner et al., 2004;Dore, 2005;Liu et al., 2013a). Here we define arid drylands as all areas that have a ratio of longterm mean annual precipitation to mean annual potential evaporation of 0.1 < P/ET p ≤ 0.3 and semi-arid drylands as 0.3 < P/ET p ≤ 0.7. Globally the primary drivers of vegetation dynamics in drylands include: (i) climate (Herrmann et al., 2005;Bai et al., 2008); (ii) fire regime (Bond and Keeley, fire, grazing and agriculture can be considered the main disturbance processes that affect large-scale vegetation trends in drylands. Other factors, e.g., nitrogen deposition and changes in growing season length (Bai et al., 2008) may also play an important role locally. Changes in the five primary drivers are expected to result in differential responses of woody cover (i.e., woody encroachment or decline) and herbaceous cover (Archer et al., 1995;Van Auken, 2000;Bond et al., 2003;Asner et al., 2004).
Although supporting evidence has been found that each of the five primary drivers can result in responses of dryland vegetation cover, simultaneous changes of the five primary drivers and interactions among them and interaction with ecosystem processes make it difficult to attribute change to any single driver (Bond and Midgley, 2012). Using remote sensing data the evidence of climate impacting vegetation dynamics has been studied at regional (Evans and Geerken, 2004;Herrmann et al., 2005;Wessels et al., 2007), continental (Donohue et al., 2009) and global scales (Nemani et al., 2003). To date, the impact of grazing and fire have not been included in global studies of dryland vegetation dynamics (e.g., Nemani et al., 2003;de Jong et al., 2011). Remote sensing data on global fire regimes are available, with many improved products becoming available since the launch of the first Moderate Resolution Imaging Spectroradiometer (MODIS) in 1999 (Kaufman et al., 1998). Studies of fire-induced vegetation change are mainly performed regionally (e.g., Flannigan et al., 2009;Archibald et al., 2010;Heisler et al., 2003) with notable exceptions focusing on longer timescales (Bowman et al., 2009;. Vegetation indices have been widely used as indicators for advances in agricultural practice (e.g., irrigation and fertilization) and annual crop yields (e.g., Tottrup and Rasmussen, 2004), while land cover maps provide information on the spatial extent of agriculture. Several studies provide field evidence of the impact of grazing (Van Auken, 2000;Asner et al., 2003) and rising CO 2 levels (McMahon et al., 2010;Buitenwerf et al., 2011) on vegetation cover; however, these studies are regional, presumably to avoid problems of complexity and data availability. As a consequence, at a global scale the relative importance of climate, fire, grazing, agriculture and CO 2 on vegetation dynamics is still actively debated and so far unresolved (Archer et al., 1995;Asner et al., 2004;Bond and Keeley, 2005;Fensham et al., 2005;Sankaran et al., 2005;Buitenwerf et al., 2011;Lehmann et al., 2011;Bond and Midgley, 2012).
Satellite remote sensing is a powerful tool to globally study both woody encroachment and desertification over multiple decades, complementing local to regional field evidence and process studies. For global studies, the Normalized Difference Vegetation Index (NDVI) is the most widely used spectral vegetation index and is based on a ratio of red and near-infrared reflectance (Rouse et al., 1974;Tucker, 1979;Beck et al., 2011). NDVI has been used as an indicator of vegetation productivity (Tucker et al., 1985;McVicar and Jupp, 1998); and is related to leaf area index (LAI; Wang et al., 2005), canopy cover fraction and the fraction of absorbed photosynthetically active radiation (fPAR; Asrar et al., 1984;Carlson and Ripley, 1997;Lu et al., 2003). Complementary to the NDVI, a recently developed data set is Vegetation Optical Depth (VOD). VOD describes the transparency of vegetation in the microwave domain and is mostly sensitive to vegetation water content (Kirdiashev et al., 1979). Owe et al. (2001) developed the Land Parameter Retrieval Model (LPRM) to derive VOD from low-frequency passive microwave observations. This model was further improved by Meesters et al. (2005) and has been applied to a series of passive microwave sensors (Owe et al., 2008). VOD is sensitive to vegetation water content in both the woody and leafy vegetation components, and provides a measure of aboveground biomass (Liu et al., 2011a).
Separately, both reflective and microwave products have been used to study vegetation dynamics. Pettorelli et al. (2005) review ecological studies that used NDVI to examine climate-and human-induced vegetation change (e.g., Evans and Geerken, 2004;Herrmann et al., 2005;Wessels et al., 2007), global vegetation trends (e.g., de Jong et al., 2011; and land degradation (e.g., Bai et al., 2008). Long VOD time series were only recently developed, and have been used to study vegetation phenology (Shi et al., 2008;Jones et al., 2011Jones et al., , 2012 and to show the impact of El Niño-Southern Oscillation on Australian vegetation cover (Liu et al., 2007). Global trends in VOD were shown to correspond to changes in precipitation, livestock (e.g., overgrazing), crop production, deforestation and fires (Liu et al., 2013a, b).
To date, only an introductory assessment of potential insights from a combined interpretation of NDVI and VOD cotrends has been performed (Liu et al., 2011a). While many regions had similar co-trends regional differences were illustrating that VOD provides new information for mixed woody-herbaceous land cover types. Due to the characteristics of both products, differential responses to changes in land cover are expected (Shi et al., 2008). Compared to VOD, NDVI saturates at relatively low biomass and is therefore most sensitive to vegetation covering the largest surface (Tucker, 1979). As a consequence, temporal correlation between NDVI and VOD is high for grass and croplands but lower for high biomass vegetation types where NDVI saturates (Liu et al., 2011a). Grasses are the main source of large interannual variation in NDVI for savanna ecosystems (Archibald and Scholes, 2007;Roderick et al., 1999b;Lu et al., 2003;Donohue et al., 2009). VOD has a greater penetration capacity and is more sensitive to changes in woody vegetation (Liu et al., 2011a;Shi et al., 2008). The relationship between NDVI and VOD is further explored in the background theory section. This paper uses the two complementary remote sensing data sets and analyzes their trends across the global drylands. In particular, we address the following questions: 3. What are the remaining trends in both vegetation indices, and how are they related to the other primary drivers (fire, grazing, agriculture and CO 2 )?
Differences between NDVI and VOD are explored and their trends and co-trends are interpreted ecologically over the global drylands. A model is developed to estimate the dryland vegetation responses that can be explained by precipitation variation. Residual trends (i.e., the observed minus model-explained trends) and their potential drivers are discussed and attributed where possible, by comparison with independent data sets (e.g., burned area and grazing) and information from previous studies.

Background theory
A conceptual framework is required to relate temporal patterns and trends in NDVI and VOD to vegetation characteristics. For NDVI, the reflectance observations from which it is derived are often described by a linear mixture of the constituent surface reflectance components (cf. Roberts et al., 1998): where ρ is the wavelength reflectance and f the fraction canopy cover, with subscripts denoting overstory (o), understory (u) and soil surface (s). This is a somewhat simplified approach, as light reflecting from one component can affect reflectance of another component by multiple scattering and transmittance (Roberts et al. , 1993). Homogeneous vegetation cover is assumed such that the overstory will obscure some part of the understory, and it is assumed that the NDVI observations are representative for nadir measurements such that f equates to projected canopy cover. Here it is assumed that the overstory and understory can be assumed equal to the (woody) persistent and (herbaceous) recurrent vegetation components (cf. Roderick et al., 1999b;Lu et al., 2003;Donohue et al., 2009). NDVI is calculated as a reflectance difference index and therefore, using Eq. (1) it can be shown that NDVI will respond to surface component mixing ratios in a slightly non-linear manner, although NDVI still converges to end member NDVI values when these dominate overall reflectance (Asner, 1998). NDVI per unit canopy area for recurrent vegetation tends to exceed that of persistent vegetation, with recurrent vegetation having shorter-lived and "greener" (i.e., higher chlorophyll density) leaves (Reich et al., 1997(Reich et al., , 2003. As a result, all else being equal, an increase in NDVI can be explained by an increase in total canopy cover or a relative increase in recurrent vegetation canopy cover, or both. VOD (denoted τ ) can be interpreted as being linearly related to total aboveground biomass (AGB) water content, i.e., the sum of water in woody and non-woody vegetation (Jackson et al., 1982;Wigneron et al., 1993): where θ is vegetation water content, m the AGB, and c τ the constant of proportionality. θ will normally be greater for herbaceous vegetation than for woody vegetation (Roderick et al., 1999a(Roderick et al., , 2000, therefore increases in VOD can mean an increase in total AGB, an increase of the fraction of herbaceous vegetation AGB, or both. Finally, the relationship between NDVI and VOD is influenced by the connection between AGB (m) and canopy cover (f ). This relationship can be presented by considering the commonly used light extinction equation (Monsi and Saeki, 2005) that relates f to leaf area index ( ): where κ is the extinction coefficient, α SLA the specific leaf area (i.e., area per unit leaf mass), and f leaf the fraction of leaf biomass in total AGB. While κ depends on leaf orientation and clumping, f leaf and α SLA will be greater for non-woody vegetation than for woody vegetation. Overall, total AGB remaining unchanged, the sensitivity of NDVI to a relative increase in the herbaceous AGB is expected to be considerably greater than the sensitivity of VOD. VOD, being sensitive to AGB water content, is expected to be more sensitive to a relative increase in woody vegetation. Per unit mass the tree foliage will have a somewhat higher relative water content than the woody parts, but because woody biomass typically represents 90 % or more of total aboveground biomass (Northup et al., 2005, and references therein) it still contains most of the water (Sternberg and Shoshany, 2001). By using the trend information from the different vegetation data sets the above considerations lead to the following four expectations: 1. a long-term increase in both NDVI and VOD signifies an increase in the relative fraction of herbaceous AGB and/or an increase in total AGB; 2. a long-term increase in NDVI combined with a decrease in VOD signifies an increase in the relative fraction of herbaceous AGB; 3. a long-term decrease in NDVI and an increase in VOD signifies an increase in the relative fraction of woody AGB; and 4. a long-term decrease in both NDVI and VOD signifies a decrease in the relative fraction of herbaceous AGB and/or a decrease in total AGB. www.biogeosciences.net/10/6657/2013/ Biogeosciences, 10, 6657-6676, 2013 To study global long-term dryland vegetation dynamics, both NDVI and VOD data sets were used for their overlapping period 1988-2008. Both vegetation indices are dimensionless and here we use a 0.25 • spatial and monthly temporal resolution.

Normalized Difference Vegetation Index (NDVI)
The NDVI data are Global Inventory Modeling and Mapping Studies 3rd generation (GIMMS3g) derived from the Advanced Very High Resolution Radiometer (AVHRR) sensor on board the NOAA series of satellites (i.e., NOAA 7, 9, 11, 14, 16 and 17;Tucker et al., 2005). Although datasets are merged from several sensors and corrections were performed , they do not affect the calculated trends (Fensholt and Proud, 2012). These data were resampled from 1/12 • to 0.25 • spatial resolution and from twice to once monthly temporal resolution by simple averaging. Data were available between 1981-2010. Recently Beck et al. (2011) showed that of four AVHRR processing chains tested, GIMMS was best able to track trends in 1424 Landsat image pairs.  tem, 2002-2008) using the LPRM algorithm. This harmonized data set preserves the relative dynamics (e.g., seasonality and interannual variations) of the original products (Liu et al., 2011b) and is able to capture long-term changes in total aboveground vegetation water content and biomass over various land cover types globally (Liu et al., 2013a). Because of enhanced sensor characteristics (particularly the longer wavelength), the accuracy of VOD retrievals from AMSR-E is expected to be better than those from TMI and SSM/I. A comprehensive evaluation study by Liu et al. (2013a) demonstrated that the errors associated to sensor changes in the harmonized VOD data set are small, however. The harmonized data set captured long-term changes in total aboveground vegetation water content and biomass over different land cover types without sensor artefacts. The VOD data should be interpreted with caution over sandy deserts, open water and under the frost conditions (de Jeu, 2003;Jones et al., 2011;Gouweleeuw et al., 2012).

Precipitation
Precipitation data (monthly; 1901-2009) produced by the University of East Anglia Climatic Research Unit (CRU Time Series version 3.1) were available at a 0.5 • resolution (Jones and Harris, 2008). For analysis it was assumed that within this grid cell precipitation was homogeneously distributed and data was resampled to a 0.25 • resolution (Koenig , 2002).

Land cover
The MODIS Land Cover Type Yearly Climate Modeling Grid (MCD12C1; 2005; 0.05 • resolution), land cover map was used here (NASA, 2008). Data was resampled to a 0.25 • resolution using the land cover with the highest frequency in each grid cell. The land cover types follow the University of Maryland (UMD) classification scheme (Hansen et al., 2000). For our analysis some land cover classes were merged to assist interpretation (i.e., urban and built-up, barren or sparsely vegetated, and unclassified + fill were here merged into "barren or sparsely vegetated"; all types of forest into "forest"; and closed and open shrublands into "shrublands").

Fire
Monthly burned area data was available based on MODIS Terra satellite imagery MCD64A1; 500 m resolution; November 2000 onwards (Giglio et al., 2009). Data was rescaled to a monthly 0.25 • resolution by calculating the mean burned area per 0.25 • grid cell.

Relationship between NDVI and VOD
To illustrate the theoretical relations between NDVI and VOD (Sect. 2), we calculated global mean values and standard deviation for both vegetation indices. To provide more insight, we calculated annual mean values and range (between maximum and minimum) and compared seasonal patterns of three southern African land cover types with increasing woody cover: (i) grasslands, (ii) savannas and (iii) woody savannas. In addition to earlier explorations on the Biogeosciences, 10, 6657-6676, 2013 www.biogeosciences.net/10/6657/2013/ relationship between NDVI and VOD (Liu et al., 2011a;Shi et al., 2008), the results of these analysis were used to verify the four expectations of Sect. 2, regarding the co-trends between NDVI and VOD. Trends in NDVI and VOD were analyzed separately (Sects. 4.2 and 4.3), but co-trends (NDVI and VOD) were also interpreted using our theoretical framework.

Response of vegetation index anomalies to antecedent precipitation anomalies
Following Evans and Geerken (2004) and Herrmann et al. (2005), linear regression models between antecedent precipitation and each vegetation index were developed for each grid cell. These regression models were used to estimate the long-term vegetation trends expected due to precipitation patterns alone. Estimated vegetation trends were calculated in three steps. Firstly, anomalies were calculated both for precipitation and the vegetation products (NDVI and VOD). Second, the antecedent precipitation index (API) was calculated, here defined as the optimal correlating precipitation running mean (PRM) as judged from Spearman's ranked correlation coefficient (R s ) over 1-to 60-month averaging periods. The API(x, y, t) was calculated for each grid cell (x,y) and month (t) using where T is the averaging period in months for antecedent precipitation that leads to the highest R s between the PRM and the vegetation anomalies and d is the precipitation anomaly. The estimated vegetation variation was calculated as a linear function of API: estimated(x, y, t) = a 1 (x, y) · API(x, y, t) + a 2 (x, y).
For each grid cell the coefficients a 1 and a 2 were determined by least squared differences between the API and vegetation anomaly time series. Similar approaches were used by Evans and Geerken (2004), Herrmann et al. (2005) and Wessels et al. (2007) for NDVI, who also assumed (as we do) that if estimated vegetation variation based on precipitation variation was removed from observed vegetation anomalies, the residual trends are largely independent of precipitation.
Here we develop models based on precipitation and vegetation anomalies (rather than the original vegetation indices) as all trends are present in the anomalies, and not in the seasonal pattern. So a model optimized for anomalies will give a better estimation of trends in the vegetation indices caused by precipitation. Although direct correlation between precipitation and vegetation indices will be higher than correlation between vegetation index and precipitation anomalies, a model based on the original vegetation indices and precipitation would be more suitable to study the effect of precipitation on seasonal rather than interannual vegetation dynamics. However, to facilitate comparison with previous studies, the analysis was also repeated using the original vegetation indices rather than anomalies (cf. Herrmann et al., 2005).

Global dryland vegetation trends
After calculating the estimated index anomalies, trends were calculated for the observed, estimated and residual vegetation index. The residual trend was calculated as the trend in the observed minus estimated vegetation anomaly and represents the component of the signal that could not be directly attributed to precipitation variations. A conventional nonparametric Mann-Kendall trend test was used to determine areas of significant monotonic trends (cf. de Jong et al., 2011;Liu et al., 2013a). The non-parametric Theil-Sen estimator of slope is insensitive to outliers and was used to calculate linear trends (Sen, 1968;Theil, 1950).
Climate and land cover affect NDVI and VOD dynamics and will likely cause them to respond differently to the primary drivers of dryland vegetation dynamics. Hence, to stratify our results, we used global maps of land cover ( Fig. 1a) and P/ET p classes; referred to as humidity classes hereafter (Fig. 1b). Global maps of livestock density, burned area, recent trends in burned area, and ecosystem characteristics of land cover and humidity were used to interpret results. Only static livestock density and land cover data were available, but these could in some cases be combined with regional studies to interpret the influence on vegetation dynamics. It was assumed that globally consistent trends not explained by climate, fire, grazing and agricultural developments are caused by increasing atmospheric CO 2 concentrations (Bond et al., 2003;Donohue et al., 2013). Impact of fire on vegetation indices was further explored studying time series of NDVI, VOD and burned area for grassland, savanna and woody savanna in southern Africa.
Grid cells with less than 40 % valid data (i.e., 100 or less of the 252-month series) in either vegetation data set were not included in this analysis. NDVI cannot be used over snow and ice (Brown et al., 2006), while VOD is sensitive to frost conditions (see Sect. 3.2; Liu et al., 2011a), hence seasonally recurrent data gaps exist in both products during winter. The threshold of 40 % valid data was chosen to include most drylands at high latitudes/elevations. Data gaps were ignored in trend calculations yet are considered in their interpretation.

Relationship between NDVI and VOD
Globally, for the arid drylands, mean VOD values were generally lower than NDVI, with increasing humidity and biomass VOD increased faster than NDVI ( Fig. 2a and b).
Most notable examples of large differences between NDVI and VOD were forests and agricultural regions of the temperate northern hemisphere drylands. NDVI generally showed Table 1. Percentage of global drylands with significant (p < 0.05) correlation between API and the two vegetation indices for land cover and humidity classes. The average humidity (P/ET p ) for each land cover class is shown in brackets with the land cover classes ordered in increasing humidity. 1 % represents ≈ 750 grid cells (0.25 • resolution), and abbreviations of land cover classes are explained in the legend of Fig. 1a   precipitation (P) data (Hijmans et al., 2005) and potential evaporation (ET p ) data (Zomer et al., 2008). All maps in this paper (except for   precipitation (P) data (Hijmans et al., 2005) and potential evaporation (ET p ) data (Zomer et al., 2008). All maps in this paper (except for Fig. 3d) are projected in the Miller cylindrical projection (60 • N-60 • S, 130 • W-160 • E). Terrestrial areas that are not drylands (0.1 < humidity ≤ 0.7) are from now on masked grey in all figures and are excluded from analysis, and oceans are masked (white) in all figures.
higher standard deviations than VOD, especially for drylands that face cold winters ( Fig. 2c and d). Regional exceptions were observed, with standard deviation of VOD being larger for some savanna regions in northern and southern Africa and Australia. More detail was provided by studying three land cover types of southern Africa (study areas are shown in Fig. 3d). With increasing woody cover, VOD increased faster than NDVI (Fig. 3a, b and c). The annual average VOD for grassland was 0.41 and for woody savanna this increased to 0.70 (i.e., a 0.29 increase). In contrast, the annual average NDVI only increased by 0.19 (i.e., from 0.43 for grassland to 0.62 for woody savanna). The NDVI range (annual maximum minus minimum) exceeded the VOD range for all three land cover types (Fig. 3a, b and c).

Response of vegetation index anomalies to antecedent precipitation anomalies
The strongest correlation (R s ) between vegetation anomalies and antecedent precipitation index (API) were observed over arid drylands (0.1 < P/ET p ≤ 0.3; Fig. 4a and c). While NDVI and VOD showed similar spatial patterns, VOD showed higher correlation coefficients, suggesting that VOD reacts stronger to interannual precipitation variability in drylands than NDVI. In regions with low or insignificant correlations between vegetation and API, factors other than precipitation are likely to determine interannual vegetation variability and/or variability was minimal over the study period. While in general the strongest correlation was found between VOD and API, the total area of significant correlation was similar for both vegetation indices (Table 1).  Averaging periods (denoted T in the methods; Fig. 4b and d) are related to the capacity of vegetation to use antecedent precipitation and the lead time of interannual variation in precipitation followed by a vegetation response. NDVI generally has shorter averaging periods than VOD (cf. Fig. 4b with d). To facilitate comparison with other studies, model results based on original vegetation indices and precipitation are shown in the annex material (Fig. A1). In general, observed correlation was higher, and averaging periods were shorter for the model based on original vegetation indices and precipitation (Fig. A1) than when compared to the model based on anomalies used here (Fig. 4, see Sect. 4.2). Figure 5 shows box plots of the data in Fig. 4, stratified by land cover and humidity classes. Both NDVI and VOD showed longest averaging periods for areas dominated by woody vegetation and croplands, shorter periods were observed for grasslands and savannas. Figure 6 shows an example of estimated and observed anomalies in NDVI and VOD for two 5 • × 5 • regions (Fig. 3d): (i) eastern Australia (25-30 • S, 145-150 • E); and (ii) southern Africa (25-30 • S, 20-25 • E). Modeled NDVI is based on shorter averaging periods, because the signal showed little interannual variation. VOD, on the other hand, is usually based on longer averaging periods and therefore follows the interannual variation in precipitation, being less sensitive to small intra-annual variability.

Global dryland vegetation trends
For some regions, observed annual trends in NDVI andVOD (1988-2008) showed similar trends: e.g., India, North America, eastern Australia, and northern African savannas ( Fig. 7a  and b). Other areas showed contrasting trends between the data sets: e.g., southern Africa, northern Australia, Argentina and central Asia. While observed trends in NDVI and VOD differ considerably, the estimated trends in NDVI and VOD were relatively similar and strong decreasing trends were found in southeastern Australia and Mongolia, while increasing trends were found in most of Africa and northern Australia ( Fig. 7c and d). Figure 7e and f show residual trends, calculated as observed anomalies minus estimated vegetation anomalies. For the NDVI (Fig. 7e) some trends persisted (e.g., positive trends over most of India and Spain, and negative trends over southern Russia and Kazakhstan). In other cases, if observed and estimated trends had the opposite direction (i.e., one is positive and the other negative) this resulted in enhanced residual trends (e.g., Argentina, arid northern Africa, southern Africa and northern Australia in the case of NDVI; compare Fig. 7a, c and e). If observed and estimated trends had the same direction this resulted in smaller residual NDVI trends (e.g., Mongolian steppe and semi-arid drylands of northern Africa), or trends tended to zero (e.g., southern India and northern China). Australia and southwest Western Australia). Observed and estimated trends in VOD generally had the same direction, resulting in smaller residual trends (e.g., northern Australia, Sahel, southern Africa and southern Argentina; compare Fig. 7b, d and f). Although there were also regions with trends in observed VOD, and no corresponding trends in estimated vegetation anomalies (e.g., east Africa and northern India). Figure 8a-c shows the co-trends between NDVI and VOD using the four categories of our conceptual framework, and their corresponding ecological interpretation is shown in Fig. 8d. Median observed trends in NDVI were increasing for savanna and croplands, and decreasing for the other land cover classes (Fig. 9a). For VOD decreasing median observed trends were found in forested drylands, while increasing median trends were found for all other land cover classes (Fig. 9a). Over arid drylands the median VOD trend was increasing while NDVI showed mostly decreasing trends, for more humid drylands NDVI and VOD trends were more similar and closer to zero (Fig. 9d). Savannas and shrublands showed a small increase in estimated median NDVI and VOD (Fig. 9b), indicating that part of the change in median vegetation index response is explained by trends in precipitation. There were no particular humidity classes that showed a large change in estimated median vegetation response, indicating that median changes in precipitation were close to zero for each humidity class (Fig. 9e). In the residual trends (i.e., observed minus estimated) similar but smaller median trends remained, both for NDVI and VOD ( Fig. 9c and f).
Vegetation trends (Fig. 7), are expected to be related to trends in precipitation (Fig. 10a) and trends in annual burned area (Fig. 10b). Most global drylands experienced stable or increasing precipitation amounts during 1988-2008; notable exceptions are southeast Australia, northern India, the Mongolian steppe and northern China (Fig. 10a). Trends in annual burned area declined for northern Africa and increased in most of southern Africa (2001-2011; Fig. 10b). Southern America showed a mixed pattern of increased and declined trends of annual burned area, while trends of annual burned area over most of Australia have been stable or declined. The relationship between burned area and vegetation indices was further analyzed using three study areas with different land cover in southern Africa (locations are shown in Fig. 3d). For grasslands, annual burned area is relatively low and no obvious relationship between vegetation indices and burned area is present (Fig. 11a). In savannas, secondary to precipitation, both NDVI and VOD are related to interannual variation in burned area (Fig. 11b; compare 2006 and 2008 with the other years). For woody savannas, increased annual burned area together with relatively low annual minimum NDVI values were observed in 2004 and 2007-2009, but interannual variation was much smaller.

Relationship between NDVI and VOD
Our conceptual framework resulted in four expectations corresponding to the four categories of Fig. 8d. In addition to exploring the relation between NDVI and VOD theoretically, we illustrate these expectations by comparing global distribution of mean and standard deviation of both vegetation indices. Following Liu et al. (2011a), we found that VOD showed stronger increase than NDVI, moving from low to high biomass regions (Figs. 2a, b and 3). Standard deviations of the vegetation indices were mostly influenced by seasonal fluctuations in greenness and AGB, driven by precipitation and/or temperature. While NDVI is more sensitive to seasonal greening in regions dominated by herbaceous vegetation (Archibald and Scholes, 2007;Donohue et al., 2009, Fig . 3), for southern Africa the seasonal biomass fluctuations seem to be better captured by the VOD signal. Ovington et al. (1963) show that seasonal AGB variation (mostly caused by herbaceous vegetation and tree leaves) for savanna vegetation is about three times larger than for grassland. This corresponds to a 2.2 times increase in VOD range from grasslands to savannas. By comparison, the NDVI range only increased 1.3 times ( Fig. 3a and b). As the exact physical relationships between NDVI and VOD are currently being resolved, we limit our theoretical framework (Sect. 2; Fig. 8) to contrasting trends in NDVI and VOD, comparing directions rather than the magnitudes of trends. Furthermore, the conceptual framework applies to both temporal and spatial vegetation dynamics. For example, an increase in both NDVI and VOD signifies an increase in the relative fraction of herbaceous AGB and/or an increase in the total AGB; this can be caused both by increased growing season length and/or by actual changes in vegetation density.

Response of vegetation index anomalies to antecedent precipitation anomalies
As expected, arid drylands showed the strongest correlation between the vegetation index anomalies and API (i.e., median R s NDVI > 0.3, median R s VOD > 0.5; see Fig. 5). Nemani et al. (2003) mapped these areas as being water limited, and for higher latitudes, water-and temperature limited. Correlation is less strong than found when analyzing the original indices ( Fig. A1; Herrmann et al., 2005) because seasonal vegetation-precipitation responses are not included. In semiarid drylands, regions that have a strong seasonal precipitation response do not necessarily show a strong interannual response (compare Fig. 4a and c with Fig. A1). This can be explained by a seasonal abundance of water in which variation in precipitation does not affect vegetation, followed by a dry season in which the vegetation is unable to use antecedent precipitation.
Following previous studies (Evans and Geerken, 2004;Herrmann et al., 2005;Nemani et al., 2003;Wessels et al., 2007;Donohue et al., 2009), we also found that not all interannual variation can be explained by precipitation alone. Precipitation-driven estimated trends in NDVI and VOD largely agree with trends in soil moisture reported by Dorigo et al. (2012). VOD generally shows longer averaging periods than NDVI, this is because they are sensitive to different components of vegetation cover (Fig. 3a-c); NDVI is more sensitive to changes in the shallow rooted herbaceous understory (Archibald and Scholes, 2007), and VOD more sensitive to changes in the woody overstory, which can also utilize moisture from deeper soil and groundwater stores (Rossatto et al., 2012;House et al., 2003). Clear geographical patterns are present in the NDVI and VOD averaging periods (see Fig. 4b and d), suggesting that physical landscape (including landform and soil type), affecting water holding capacity, likely contributes. The stronger correlations between VOD and API are attributed to the larger interannual variation in the VOD data set and caused by interannual precipitation variations (Fig. 6). Archibald and Scholes (2007) reported similar findings, and concluded that plants that access deeper water and have carbohydrate reserves may show a phenology that is quite different from surrounding areas with grass cover that depend on shallow soil moisture for their growth. It appears that similar conclusions may be drawn for interannual variation.

Global dryland vegetation trends
The relationship between NDVI and VOD trends (Fig. 8) provides new insights in the relative performance of herbaceous and woody vegetation components in global drylands. Woody encroachment into grassland or savannas has been observed across the globe, including Argentina (e.g., Adamoli et al., 1990), Africa (e.g., Vegten, 1984;Oldeland et al., 2010) and Australia (Fensham et al., 2005, cf. review by Archer et al., 2001). Although satellite observations of NDVI and VOD confirm our theoretical framework, future validation studies, comparing satellite observations with local to regional studies, would help to improve  regional interpretation and confidence. Vegetation trends of the world's arid drylands and semi-arid drylands are discussed in turn. Overall, the world's drylands show a decrease in NDVI and an increase in VOD (Fig. 9d).

Arid drylands (0.1 < P/ET p ≤ 0.3)
The world's arid drylands are mainly covered by shrublands (48.2 %) and grasslands (31.3 %). Trends over bare or sparsely vegetated areas (10.6 %) need to be treated with care due to the limitations of both vegetation indices over those regions (see Sect. 3;Tucker et al., 2005;Brown et al., 2006). Many arid drylands experienced unchanged or increasing annual precipitation (Fig. 10b), and median estimated trends in both vegetation indices were close to zero (Fig. 9e). NDVI median trend for arid drylands was negative while VOD median trend was positive (Fig. 9e). Both the relatively constant or increasing VOD trends and most negative trends in NDVI remain after precipitation-induced variation is accounted for ( Fig. 7e and f). This results in many regions showing opposite trends in NDVI and VOD ( Fig. 8; e.g., arid drylands of Argentina, southern Africa, northern Africa and Australia).
In those regions, it seems that woody encroachment takes place at the expense of the herbaceous understory. Woody encroachment has been attributed to climate (Fensham et al., 2005), fire regime (Sankey et al., 2012), grazing (Asner et al., 2004) and CO 2 fertilization (Buitenwerf et al., 2011). Clear trends in the NDVI and VOD residuals indicate that next to precipitation, other drivers also play a role ( Fig. 7e  and f). A possible explanation for residual trends could be changing fire regimes (Fig. 10b) through its impact on competitiveness of the herbaceous and woody vegetation components (Sankey et al., 2012). Bowman et al. (2009) showed that annual burned area is highest in areas of intermediate primary production, limited by a lack of dry periods towards the tropics and limited by a lack of fuel towards dry areas (Fig. 10d). Primary production is the main limitation on annual burned area in these drier regions (Archibald et al., 2009) and increasing precipitation and/or CO 2 concentrations (Donohue et al., 2013) can increase net primary production and hence increase annual burned area in arid drylands. This results in a system of highly variable biomass production, followed by infrequent fires (Fig. 11a). Recent trends (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011) in annual burned area do not show increases in most arid drylands and showed declines in northern Australia and the Sahel (Fig. 10b). This, combined with an average annual burned area of below 10 % (Fig. 10d), seems to make fire an unlikely driver of globally observed changes.
Decreasing NDVI trends have been interpreted as a proxy for land degradation (Wessels et al., 2007;Bai et al., 2008) which might be caused by grazing. Grazing has also been associated with woody encroachment (e.g., Asner et al., 2004). Although grazing might cause contradicting trends between NDVI and VOD in some regions (e.g., arid parts of northern Africa) similar trends are observed in regions of limited or no grazing by domestic livestock (e.g., compare Figs. 10c and 8 for arid drylands of Argentina and Australia). The spatial patterns and scale at which trends occur suggest that grazing by domestic livestock is not the main driver in these cases, impact of grazing by non-domesticated animals was not included here, as no information was available. Sankaran et al. (2005) suggest that woody vegetation receiving less than 350 mm annual precipitation is largely constrained by water availability. Rising atmospheric CO 2 concentrations affects water use efficiency and photosynthetic rates (Donohue et al., 2013;Franks et al., 2013;Keenan et al., 2013) as well as light and nutrient efficiency (Drake et al., 1997;Farquhar, 1997). CO 2 is reported to enhance the relative performance of woody C3 species over the C4 grasses that dominate tropical savannas (Bond and Midgley, 2012;Higgins and Scheiter, 2012). Once woody plants are established in savannas, they are likely to limit the growth of herbaceous plants through water uptake and shading (Breshears, 2006).
Increases in extreme air temperatures and changes in growing season duration can also explain some of the observed vegetation trends (Allen et al., 2010). Changing duration of rainy seasons can impact the relative performance of different species: a shorter but more intense wet season is likely to result in decreasing herbaceous understory (NDVI) during dry seasons, while deep rooted vegetation suffers less, having access to deeper water resources and therefore being more sensitive to long-term water availability (Fig. 5;Archibald and Scholes, 2007). Evidence for shortening growing seasons is found for arid drylands of northern Africa (de Jong et al., 2011). Finally, increased air temperatures can result in higher evaporation rates and during the dry season result in a competitive advantage of deep rooted species over shallow rooted species (Tietjen et al., 2010). Although air temperature is increasing over many drylands globally, decreasing trends are also observed in some regions, as temperature tends to be strongly related to precipitation through the effect of cloud cover on solar irradiation. Additionally, decreases in wind speed ("wind stilling"), and changes in other meteorological variables governing the evaporative process, also alter the evaporative regime (McVicar et al., 2012). Following Higgins and Scheiter (2012), we suggest that the effects of changing evaporation rates in drylands are less important than the effects of increasing atmospheric CO 2 concentrations. Given that patterns in the residual trends are found in most global arid drylands, all experiencing very different fire and grazing conditions (compare Fig. 7e and f with Fig. 10c and d), a global driver such as CO 2 fertilization appears more plausible. trends in vegetation indices differs between humidity classes and land cover classes (Fig. 9). The median trend in NDVI is around zero for semi-arid drylands, showing positive trends over savannas and croplands and negative median trends for other land cover classes. VOD shows positive median trends for all land cover classes except forest (see Fig. 9a and d).

Semi
In most semi-arid areas, where savanna grasslands dominate, resources are available to support forests which are suppressed by frequent fire (Bond and Keeley, 2005). Fire is thought to be the next most important driver of vegetation variation in savannas after precipitation (Sankaran et al., 2005). Although northern and southern African semi-arid drylands experience similar increasing trends in precipitation (Fig. 10a), the changes in the vegetation indices are very different. In southern Africa (i.e., south of 5 • S), NDVI trends are declining, while VOD trends are increasing-likely caused by an increase in the relative fraction of woody AGB. Northern African savannas (in contrast with the grasslands in adjacent arid drylands), show increasing trends in NDVI, along with increasing trends (in most regions) for VOD; interpreted as an increase in the relative fraction of herbaceous AGB and/or increase in total AGB. The increasing trend in the northern African grasslands and savannas is widely discussed in literature and has been explained as a result of recovery after drought (1983)(1984)(1985) being the driest years; Anyamba and Tucker, 2005) and improved land management (e.g., irrigation, and soil and water conservation; Herrmann et al., 2005). Strong trends in NDVI occur in semi-arid areas with large annual burned area (compare Fig. 7e with Fig. 10d). Northern African savannas, with recent declines in annual burned area (Fig. 10b), showed increasing NDVI values (Fig. 7e). Widespread grazing (Fig. 10c) and associated fire suppression may explain decreasing annual burned area in northern Africa between 2001 and 2011 ( Fig. 10b; Archibald et al., 2010). Southern African regions, with recent increases in annual burned area, showed declining NDVI trends. We found that annual variation in burned area especially affects NDVI and VOD minimum values (Fig. 11). Human land use practice is also thought to be an important driver of annual burned area (Archibald et al., 2009). The effect of fire on vegetation indices is most profound in African and Australian savannas and woody savannas, where annual burned area typically ranges from 20 to 80 % (Fig. 10d). Given the high percentage of area burnt each year in many African savannas and woody savannas (Fig. 11), it seems feasible that changing fire regimes, affecting annual minimum values of NDVI, cause the observed trend in NDVI.
Outside the African savannas and woody savannas, the residual trends in NDVI are not easily attributed to changes  in annual burned area, at least not using the data available to us. While even in the less fire prone regions, fire might play a crucial role in species competition, the direct effect of the smaller relative area of fire scars on the vegetation indices would be limited. NDVI can also be affected by grazing; while any given area might support livestock during most times, the highest pressure on vegetation occurs during longer dry seasons. The negative trends in NDVI and VOD over the South American dry forests are likely caused by deforestation, that is reducing the woody vegetation component and therefore total AGB (Grau et al., 2005). Trends in VOD occurred in areas of frequent fire, but no coherent spatial pattern was apparent (compare Fig. 7f with Fig. 10b and 10d). Although fire suppresses woody encroachment, dryland fires generally have a relatively low temperature and flame height and therefore do not necessarily affect established woody vegetation (Bond and Keeley, 2005). Woody vegetation increased in the southern African savannas despite recently increased burning (Figs. 10b and 8). This may be explained by increasing precipitation (Fig. 7; Sankaran et al., 2005), CO 2 fertilization (Buitenwerf et al., 2011;Higgins and Scheiter, 2012) and/or grazing by nondomesticated herbivores. A clear increasing NDVI and VOD trend was found for croplands and adjacent grasslands in the USA. Woody encroachment in US grasslands has been widely reported (see Archer et al., 2001, and the references therein) and although there is no straightforward single driver, grazing is understood to play an important role (Van Auken, 2000;Briggs et al., 2005). Increasing agricultural activity could also play a role in explaining those trends (Neigh et al., 2008).
The strongest positive median NDVI and VOD trends were found in agricultural areas (Fig. 9). Increases of both indices over the world's agricultural regions are explained by advances in agricultural practice including mechanization, irrigation and fertilization (Liu et al., 2013a). In India, Pakistan, Bangladesh, China, Ukraine, southwestern Russia, several European countries, the USA, and a number of other countries, substantial areas of agricultural land are irrigated (Wada et al., 2010). Evidence that increasing trends in NDVI are caused by irrigation and fertilization has been documented for India (Jeyaseelan et al., 2007) and the North China Plain (Piao et al., 2003), while Liu et al. (2013a) showed that positive VOD trends in southern Russia, China, India and the US are the result of increased agricultural production. We used a land cover map of 2005 and for many regions agricultural extent increased during the study period. Therefore, trends may partly be explained by land cover conversion from natural land cover classes to cropland. In drylands at higher latitudes both water and temperature are factors limiting vegetation growth (Nemani et al., 2003), and global changes in temperature can therefore also explain some of the vegetation trends (Tucker et al., 2001).

Conclusions
A recently developed passive microwave vegetation data set (VOD) and a widely used reflective based vegetation data set (NDVI) were combined to study the long-term  vegetation changes over the world's drylands. We draw the following conclusions: 1. The two data sets provide complementary information on vegetation dynamics; NDVI being most responsive to canopy cover and greenness, and VOD to aboveground biomass.
2. NDVI was more sensitive to herbaceous vegetation changes and short-term precipitation variations. VOD, on the other hand, was more sensitive to changes in woody vegetation and longer-term precipitation variations.
3. Although precipitation is an important driver for dryland vegetation dynamics, precipitation variations could not explain all of the observed trends in vegetation indices.
4. Co-trends between NDVI and VOD provide evidence of widespread woody vegetation encroachment at the expense of the herbaceous vegetation component in arid regions (humidity < 0.3), and arid shrublands in particular. Spatial distribution of trends suggests that a global driver (e.g., CO 2 fertilization) is causing a change in relative performance of woody vegetation compared to herbaceous vegetation.
5. Remote sensing evidence for woody thickening and encroachment is also found for some semi-arid drylands, but regional trends vary widely. It is interpreted that local rather than global drivers are responsible for most of the observed residual trends in these areas. Limited observations of monthly burned area suggests that after precipitation, changing fire regimes are an important driver of vegetation change in semi-arid drylands, especially in savannas.
6. Large changes in vegetation density were observed in agricultural dryland regions, where advances in agricultural practices caused increasing trends in both vegetation indices.
In summary, we demonstrated that using two complementary vegetation indices provide new insights into the dynamics of different vegetation components in global drylands. While it remains challenging to conclusively attribute dryland vegetation dynamics to any individual driver, a linear precipitation response model showed that change cannot be attributed to precipitation alone. Global data on fire regimes and grazing enabled a first assessment of the likely relative importance of these drivers on global vegetation change. Future improvements and extensions to time series of fire characteristics, grazing and land use change are likely to further improve understanding of global vegetation changes and its drivers.