Global patterns of NDVI-indicated vegetation extremes and their sensitivity to climate extremes

Extremes in climate have significant impacts on ecosystems and are expected to increase under future climate change. Extremes in vegetation could capture such impacts and indicate the vulnerability of ecosystems, but currently have not received a global long-term assessment. In this study, a robust method has been developed to detect significant extremes (low values) in biweekly time series of global normalized difference vegetation index (NDVI) from 1982 to 2006 and thus to acquire a global pattern of vegetation extreme frequency. This pattern coincides with vegetation vulnerability patterns suggested by earlier studies using different methods over different time spans, indicating a consistent mechanism of regulation. Vegetation extremes were found to aggregate in Amazonia and in the semi-arid and semi-humid regions in low and middle latitudes, while they seldom occurred in high latitudes. Among the environmental variables studied, extreme low precipitation has the highest slope against extreme vegetation. For the eight biomes analyzed, these slopes are highest in temperate broadleaf forest and temperate grassland, suggesting a higher sensitivity in these environments. The results presented here contradict the hypothesis that vegetation in water-limited semi-arid and semi-humid regions might be adapted to drought and suggest that vegetation in these regions (especially temperate broadleaf forest and temperate grassland) is highly prone to vegetation extreme events under more severe precipitation extremes. It is also suggested here that more attention be paid to precipitation-induced vegetation changes than to temperature-induced events.


Introduction
Accumulating evidence has revealed the impact of climate extremes on tree mortality (Allen et al 2010, Nakagawa et al 2000, gross primary productivity (Ciais et al 2005), forest biomass (Phillips et al 2009), and other important features of ecosystems (García-Herrera et al 2010, Xu et al Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. , Schwalm et al 2012. With the projected increase in such climate extremes under future climate change (Luber and McGeehin 2008, O'Gorman and Schneider 2009, Parry et al 2007, knowledge of the severity and patterns of such ecological impacts becomes crucial because extremes in terms of climate variables are not necessarily associated with their impact (Seneviratne et al 2012). The inadequate growing conditions that vegetation reveals could be regarded as vegetation extremes, corresponding to climate extremes, and could provide us with insights into this knowledge, revealing the underlying patterns of vegetation vulnerability and demonstrating vegetation's sensitivity to climate extremes. However, while climate extremes have been widely studied at both national and global scales , Meehl and Tebaldi 2004, Parry et al 2007, Soulé 2006, global assessments of vegetation extremes remain scarce. This may be partially due to the difficulty in formulating a uniform definition of vegetation extremes that enables both temporal and spatial comparisons.
Although understanding the vulnerability of an ecosystem is crucial to predicting its dynamics, the term 'ecosystem vulnerability' has been defined differently in different studies. Some studies have defined it as sensitivity to a given adverse event (Christensen et al 2004, Propastin et al 2010, while others have used it to indicate the likelihood of abnormal vegetation dynamics (Villers-Ruiz andTrejo-Vázquez 1997, Pressey andTaffs 2001). Revealing vulnerability through frequency of vegetation extremes is equivalent to the latter approach and will be more meaningful in comparisons with various approaches used in earlier studies. Looking at response time scales, for example, Vicente-Serrano et al (2013) evaluated the sensitivity of different biomes to drought and concluded that semi-arid and semi-humid areas might be more resilient because they tend to have longer response times. However, slow response might alternatively be the consequence of a more lasting or even profound impact that becomes stronger after a certain time lag. Therefore, one could propose the converse hypothesis that semi-arid and semi-humid areas are more vulnerable than others and test this hypothesis using the cumulated frequency of vegetation extremes over a certain period.
Remote sensing data with continuous time series on different spatial scales provide a consistent and repeatable measurement of vegetation condition that is appropriate for capturing the effects of many processes that cause changes (Verbesselt et al 2010). The NDVI (normalized difference vegetation index) is a widely used index of vegetation growing conditions (Linthicum et al 1999, Lotsch et al 2005, and its low-value events exceeding certain temporal and spatial extents have been used as a proxy for vegetation extremes in this study. Instead of focusing on the response of vegetation after climate extremes, this research has analyzed patterns of vegetation extreme frequency and their relationship with climate extreme frequency. The objectives of this research are to: (1) evaluate global vegetation vulnerability through the frequency of vegetation extremes at a spatial resolution of 0.5 • , (2) analyze the controlling climate factors of this global vegetation vulnerability pattern, and (3) determine the dominant effects of climate extremes on vegetation extremes at a global scale and comparing the relationship between climate extremes and vegetation in different biomes.

NDVI data and preliminary processing
Corrected from the Pathfinder NDVI of Advanced Very High Resolution Radiometer (AVHRR) data, the GIMMS-NDVI data set with a resolution of 8 km × 8 km has become one of the most commonly used NDVI data sets . Biweekly NDVI data from January 1982 to December 2006 were used to represent vegetation conditions in this study at a resampled resolution of 0.5 • . Snow-cover data derived from MODIS (MOD10C1, Hall and Salomonson 2000) were used to mask snow-covered areas that might present additional land-cover dynamics other than vegetation. Because snow-cover data are available since February 2000, the data from February 2000 to April 2012 were averaged, and the 75th percentile of all data (multi-year average values for all seasons and all pixels) was defined as the threshold for a snow-covered two-week period for each pixel.
The NDVI data were first processed as follows: (1) Pixels with constant negative values were excluded because the NDVI of vegetation cover should be greater than zero, while the calculated NDVI can range from −1 to 1.
(2) Any remaining pixels with negative values were reset to zero for the corresponding episode. Although a negative NDVI does not include vegetation information, its value could influence the detection of vegetation extremes.
(3) For each pixel, two-week periods with snow cover greater than the threshold described above were excluded, assuming that no vegetation extremes occurred in snow-covered periods. Although using also NDVI during the growing season could greatly reduce the influence of non-vegetation land cover, partial information on vegetation would be lost, and the variance of phenology on different sites would make vegetation extreme frequencies less comparable. Therefore, only the influence of snow, which is the most significant interference with NDVI, was removed. (4) The NDVI was resampled to a resolution of 0.5 • × 0.5 • , using the average of the highest 25% of pixels within the grid. Re-compositing NDVI data to a 0.5 • ×0.5 • grid from the original 8 km×8 km grid using maximum values could reduce atmospheric influences that nearly always decrease NDVI, but will also diminish the representativeness of the NDVI value and increase sampling uncertainty. As a compromise, the mean value of the highest 25% pixels within every 0.5 • × 0.5 • grid was used. The number of values used in re-compositing the NDVI (the highest 25% of pixels in this study) is empirical, but the uncertainty of this method was evaluated using the coefficient of variation (CV) of different re-composition methods: the coefficient of variation (CV) was calculated for the mean values of the highest 15%, 20%, 25% (used in this study), 30%, and 35% pixels in every grid cell.
All the work described below was based on these processed 25-year 0.5 • biweekly NDVI data.

Definition and extraction of extremes
The most widely used approach for estimating extreme values might well be extreme value theory (EVT), which 'aims at deriving a probability distribution of events from the tail of a probability distribution' (Seneviratne et al 2012). By predefining events that occur below a certain probability as extreme, this method identifies the same number of extremes in any distribution and is therefore inappropriate for this study. Although drought duration and seasonality are crucial to ecological impact in studies of drought events (e.g., Asner and Alencar 2010), their effect is less important in this study because the vegetation extremes are already ecological consequences of adverse effects (including drought).
The most common method of vegetation extreme extraction in previous studies is the calculation of NDVI anomalies: departures from the average normalized value by one standard deviation or more (Linthicum et al 1999, Liu and Juárez 2001, Lotsch et al 2005, Xu et al 2012. Vegetation extreme events could be detected and compared if the NDVI at all sites had the same distribution; however, this assumption is violated by the nonlinear nature of the NDVI and its saturation in high-vegetation-cover areas (Huete et al 1997, Mutanga andSkidmore 2004), making comparison between different sites and different seasons problematic.
Acquiring a global pattern of vegetation extreme frequency requires a robust statistical method that can produce comparable results on a global scale. Breaks for additive seasonal and trend (BFAST), a robust approach that integrates the decomposition of time series into trend, seasonal, and remainder components, has been developed recently (Verbesselt et al 2010). However, the vegetation extremes reflected in its remainder component have a statistical vagueness similar to that of anomalies. Independent component analysis (ICA) on NDVI time series, a relatively novel approach to capture variability arising from independent physical sources, and potentially from vegetation extreme events in a given region, has been widely used in recent years (Lotsch et al 2003, Philippon et al 2007, Antico 2012. Similarly to principal component analysis (PCA) (Lasaponara 2006), wavelet decomposition (Anyamba and Eastman 1996), and change vector analysis (CVA) (Lambin and Strahler 1994) used in previous studies, this method is based on multiple time series and requires a multi-temporal spectral space, such as the areas of vegetation extremes in this study, to be defined beforehand. Therefore, this method is also inappropriate for the current work.
A new method has therefore been developed to extract vegetation extreme events from NDVI time series. It is assumed that: (1) Vegetation type remained constant through the study period of 25 years . Formal studies have indicated that change in vegetation type occurs at a time scale of decades to centuries and can therefore be ignored.
(2) Vegetation growth has a period of 1 year. This is clearly valid for most regions and can also be accepted for vegetation controlled mainly by non-seasonal precipitation.
(3) For a given pixel and a given two-week period, NDVI values for different years have a unimodal probability distribution. For a given vegetation type, the NDVI should fluctuate around a certain value at the same time of each year under normal vegetation conditions, displaying a unimodal probability distribution.
NDVI values over the 25 years for a given pixel and a given two-week period constitute a sample set. The differences in skewness and kurtosis of the NDVI distribution make comparison with earlier studies problematic. However, with a unimodal probability distribution, data can be transformed to a normal distribution. In the present study, the Box-Cox transformation (Box and Cox 1964) was used, which transforms a non-normally distributed data set X into a normally distributed data set X through: where λ is the value corresponding to the maximum L of the log likelihood function: where v and n are the degrees of freedom and the sample size, x i is an individual in X, and s 2 is the variance of X . Each group X can be processed as described below: (I) Assuming the 25 individuals of group X (figure 1(A)) to be normal individuals, λ is calculated using (2), and X is transformed into X using (1) (figure 1(B)). (II) The farthest individual from the average is regarded as a suspect individual x s , while the rest are assumed to be normal individuals and are used to calculate λ using (2), following which the whole group (including x s ) is transformed using λ and (1) (figure 1(C)). (III) Assuming a normal distribution of the transformed data set (except for x s ), according to the relationship between Student's t-distribution and normal distribution, x s should be regarded as an extreme if: where µ and s is the mean and standard deviation of the rest of the dataset, p is the confidence level (p = 0.9 in this study) and t (p,n−2) is the pth quantile of a t-distribution with degrees of freedom n − 1. If x s is not an extreme, then all individuals are regarded as normal; otherwise x s is recorded, and steps I and II are repeated with x s excluded. For the results of step I, the Lilliefors test (Lilliefors 1967), a normalized test suitable for small sample sizes, was used to verify the effect of the Box-Cox transformation. With significance level α = 0.05, 97% of the groups passed the test. The confidence level in step III is actually lower than p = 0.9 because the estimation of λ increased the uncertainty; however, p must be relatively low to provide more extremes for later analysis. After the extremes were extracted from each multi-year biweekly group, the binary data indicating whether extremes had occurred were put back into order (600 two-week periods over 25 years).  (3)) directly would identify the yellow circle on the right as an extreme, while this individual might be normal in such a distribution. (B) Assuming all individuals to be normal individuals and applying a normal transformation to the data set, the yellow circle on the left is now most likely to be an extreme. (C) Assuming the yellow circle on the left to be an extreme, the rest of the data set was used to calculate λ , and a transformation of the whole data set before the t-test was used to decide whether the yellow circle on the left is an extreme.
In the longitude-latitude-temporal three-dimensional space, the spatial and temporal extent of each set of continuous (continuous in either space or time) extremes was then calculated, and sets of extremes that were shorter than six weeks or smaller than 1 • in longitude or latitude were eliminated. The remaining events were vegetation extremes that met the research requirements. The total number of all vegetation events under all conditions in each pixel was 600 (600 two-week periods over 25 years). The total number of vegetation extreme events in the pixel was then added up to yield a vegetation extreme frequency representing the total time length of vegetation extremes in every pixel over the 25 years.

Climate variables and relationship with vegetation extremes
Quantifying climate extremes such as drought is difficult, and therefore in this study, extremes in precipitation, temperature, and drought index were evaluated separately. Global temperature and precipitation values from a CRU3.0 raster data set with a spatial resolution of 0.5 • and a temporal resolution of one month (Climatic Research Unit, www.cru. uea.ac.uk, Mitchell and Jones 2005), as well as the Palmer Drought Severity Index (PDSI) with a spatial resolution of 2.5 • × 2.5 • and a temporal resolution of one month (Dai et al 2004), were used. Extreme events in temperature (positive extremes), precipitation (negative extremes), and PDSI (negative extremes) were extracted using the identical method. Because the temporal and spatial resolution of these measures was coarser than NDVI, spatial and temporal filtering (eliminating outliers shorter than one-and-one-half months on the temporal scale or smaller than 1 • in longitude or latitude on the spatial scale) was not used on PDSI extremes, and only the spatial filter was used on temperature and precipitation extremes.
To ensure equal land areas and biome distributions in the following analysis, all original geographic reference data were transformed into an equal-area 'sinusoidal' projection using the ArcGIS software. Pixels were grouped into climate bins using 100 equally spaced intervals of global mean annual temperature (MAT), mean annual precipitation (MAP, logarithmically transformed), and mean annual PDSI. Then the vegetation extreme frequency of each bin was calculated by averaging the extreme frequency of all pixels in the bin. For comparison with climate types in the Köppen-Geiger climate classification (Peel et al 2007), the predominant climate type within each bin was also calculated.
As vegetation vulnerability was defined as 'the likelihood of abnormal vegetation dynamics', higher climate extreme frequency tends to result in higher vegetation extreme frequency and thus higher vegetation vulnerability. However, this effect is different in different biomes, depending on the sensitivity of vegetation, and the slope of linear regression models was used to evaluate this sensitivity. Eight biomes were then reclassified based on the 14 major habitat types identified by the WWF in The Global 200: Priority Ecoregions for Global Conservation (Olson and Dinerstein 2002; table 1 and figure 2(c)). For each biome, linear regression models were computed using all pixels as samples, the climate extreme frequency (precipitation, temperature, and PDSI) as the independent variable, and the vegetation extreme frequency as the dependent variable. The slopes of these regression models in different biomes were then compared. Despite the uncertainties in the relationship between extreme events in climate and vegetation, for instance nonlinearities and time lags, pixels that accumulated more climate extremes should have experienced more vegetation extremes, and therefore the slopes of the linear regressions could be used as an indicator of apparent sensitivity. Besides linear regression using the method of least squares, quantile  Tropical and subtropical grassland Tropical and subtropical grasslands, savannas, and shrublands 8 Tropical and subtropical forest Tropical and subtropical moist broadleaf forests Tropical and subtropical dry broadleaf forests Tropical and subtropical coniferous forests Mangroves regression with τ = 0.9, 0.95, and 0.99 was also used. Unlike the least-squares method, which estimates the conditional mean of the dependent variables, quantile regression estimates the given quantiles (τ ) of the response variable (Koenker 2005). Because microclimates which were ignored at the spatial scale of this research could function as refuges under climate change (Foster 2001, Murphy andWeiss 1992) and reduce the sensitivity displayed by some pixels, quantile regressions focusing on the higher values might provide a better evaluation of sensitivity and were therefore used as well.

Global patterns of vegetation extreme frequency
The coefficients of variation of the different sampling methods described in section 2.1 were plotted to illustrate the uncertainty of the sampling methods in different regions. Apart from mountainous and coastal regions with high topographical heterogeneity, most regions showed relatively low CV values ( figure 2(B)). Regions of markedly high vegetation extreme frequency were generally distributed in the semi-arid and semi-humid regions and the Amazon rainforest, including savanna in southern Africa and the Sahel, forest-grassland ecotones in North America, Central America, South America, and southeastern Australia, and Mediterranean vegetation in southern Europe. Regions with markedly low frequencies (or distinctly lower than their surroundings) include desert (the Sahara and Namib Deserts in Africa, the Arabian Peninsula, Aral Karakum in Asia, deserts in central and northern Australia) and mountainous regions (the Rockies in North America, the Andes in South America, and East African mountains along the East African Rift) (figure 2(A)). Among the eight biomes, tundra, taiga, and desert have low average vegetation extreme frequencies, while the average vegetation extreme frequency of the other five biomes varies little (figure 2(C)).

Climate variables and their relationship with vegetation extremes
The average vegetation extreme frequencies showed no clear pattern with mean annual PDSI, and therefore only the results for MAT and MAP were plotted, yielding an average vegetation extreme frequency for each of the 100 × 100 bins (figure 3). Pixels of high vegetation extreme frequency are dense for MAT of 10-28 • C and MAP of 300-1000 mm, including all or parts of Csa (dry-summer subtropical), Csb (Mediterranean), Dfb (warm-summer continental), Cfc (oceanic), BWh (low-latitude desert), BSk (middle-latitude steppe), Cwb (dry-winter temperate highland tropical), and Aw (savanna) in the Köppen-Geiger climate classification. Regions of low MAT or low MAP have low vegetation extreme frequencies corresponding to desert, boreal forest, and tundra (figure 2). In regression analysis, each dot corresponds to one pixel, while the independent variable is the precipitation extreme frequency and the dependent variable is the NDVI extreme (vegetation extreme) frequency (figure 4). The numbers shown are the slopes of different linear regression models. Significant positive correlations (p < 0.001) were found for most biomes and most models. The slope increased with τ in quantile regression for all biomes except tropical and subtropical forest. For all four models, the slopes were highest in temperate broadleaf forest and temperate grassland. The distribution of vegetation extreme frequency in each biome (shown by the gray dots) also indicates that the highest vegetation extreme frequencies were found in temperate broadleaf forest and temperate grassland.

Spatial consistency between high frequencies of vegetation extremes and climate extremes
A global pattern of vegetation extreme frequency was successfully estimated in this study by developing a robust method of detecting significant extremes (low values) in biweekly time series of global normalized difference vegetation index (NDVI) values from 1982 to 2006. Unlike earlier studies, this method considers low-NDVI events and their spatial and temporal scales simultaneously and guarantees the same statistical significance for extremes extracted from different regions in different seasons, making vegetation extremes extracted from different sites comparable. Due to the Box-Cox transformation and Student's tdistribution, the threshold of vegetation extreme is different for different sample groups, depending on their own variation and distribution; therefore, vegetation extremes are equally unlikely events for the given pixel and biweek, and are not more frequent in places with higher NDVI variation.
Direct global-scale validation of the consequential extreme frequency pattern cannot be achieved because the method used in this study extracted seasonal vegetation extremes, unlike the traditional method used for drought events that considers duration and spatial extent (Sheffield et al 2009). However, although vegetation extremes in this study could occur more than once each year, their frequencies partially indicate a cumulative duration of vegetation extremes and are therefore comparable with previous research focusing on single extreme events, normally of indefinite duration. The vegetation extreme frequency pattern agrees well with previous regional studies using various methods, especially in Europe, America, and Africa, which have been widely studied. At a global scale, most of the forest mortality areas are situated in the regions of high frequency in the present results (Allen et al 2010), which have been largely attributed to drought and heat stress. Therefore, the results obtained here were compared mostly with severe drought and heat-wave events.
In Europe, the highest vegetation extreme frequencies appeared in France, Spain, northern Italy, Romania, Serbia, Greece, and on the boundary between Ukraine and Russia, covering Mediterranean vegetation and montane forests, or alternatively, BSk, Csa, Csb, and Cfb in the Köppen-Geiger climate classification. This pattern agrees with previous studies of both vegetation vulnerability patterns and single climate-induced vegetation decline events in Europe. Spain, Italy, and Romania are also the highest in Europe in term of vulnerability of their major agricultural areas (Gol'tsberg andPokrovskaya 1972, cited by Kogan 1997). During the 2003 heat wave in Europe, the largest estimated net primary productivity (NPP) reductions were found in Ukraine, Romania, France, and Italy (Ciais et al 2005). A similar pattern was also found in an overview of forest decline in southern Europe (Bussotti and Ferretti 1998).
In North America, vegetation extreme frequency is highest in the interior and western portions of North America, except for the Rocky Mountains and the Sierra Madre. The most severely stressed regions during the dry spells in 1988, 1989, and 1996, as indicated by the AVHRR-based NDVI anomaly and vegetation condition index (VCI), were Texas, Arkansas, Mississippi, Alabama, Ohio, and West Virginia in 1988, the southwestern part of the United States except for parts of Utah, Colorado, Arizona, and New Mexico in 1989, and in addition the southwestern region of the Great Lakes in 1996 (Kogan 1995(Kogan , 1997. Patterns of severe vegetation stress were different in these three events, but their overlap coincides with the high vegetation extreme frequency regions of North America in the result of this study. The frequency of months under drought and wet spells from 1895 to 1981 in every state, although at a coarser resolution, also exhibited a pattern similar to these results, highlighting the interior and western portions of the United States (Diaz 1983).
The high vegetation extreme frequency regions in South America include the Amazon rainforest and tropical and subtropical dry-forest regions in southern Brazil, Paraguay, and northern Argentina. With the most severely affected areas in the eastern part of Amazonia and the grassy savanna in northern Brazil, this pattern roughly coincides with regions of NPP decline from 2000 to 2009 in South America (Zhao and Running 2010). Most parts of this region have been identified as 'bistable forest' or 'bistable low tree cover', which is capable of transforming between forest and savanna (Staver et al 2011). Frequent vegetation extremes could be partially attributed to the strong influences of the tropical Atlantic north-south SST gradient on dry-season rainfall, because its intensification (warming of northern SSTs relative to the south) shifts the Inter-tropical Convergence Zone northwards (on inter-annual time scales) and intensifies the dry season in southern and eastern Amazonia (Li et al 2006), as occurred in the 2005 drought (Malhi et al 2008, Marengo et al 2008. Agricultural colonization might also have caused forest decline in this region, but was restricted to local scales, for instance the expansion of cattle and soybean production in Amazonia (Malhi et al 2008) and logging in eastern Paraguay (Riezebos and Loerts 1998).
The high vegetation extreme frequency regions in Africa are the Sahel and the tree savanna/temperate grassland areas of southern Africa. A similar pattern was acquired using the variable infiltration capacity (VIC) model to identify the 1983, 1984, and 1992 drought events in Africa (Sheffield et al 2009). The Sahel is the most dramatic example of climate variability that has been directly measured (Hulme 2001); its drought-induced vegetation declines in recent decades have triggered an upsurge in interest in related scientific research as well as political debate (Herrmann et al 2005). It has also been identified as 'bistable low tree cover' which can change to forest (Staver et al 2011). The areas with the highest vegetation extreme frequencies in the Sahel region in the current results were the eastern and western ends, which is similar to the vegetation decline pattern in the Sahel during the 1983 and 1984 droughts across Africa . The tree savanna and temperate grassland in southern Africa, including the northern and eastern parts of the Republic of South Africa, Mozambique, southern Malawi, Zimbabwe, northern Botswana, and southern Angola, were found to be among the most severely affected regions in the 1989, 1992, and 1993 droughts in southern Africa in terms of the vegetation condition index (VCI) and the temperature condition index (TCI) (Unganai and Kogan 1998).
These areas of agreement with previous studies indicate that the proposed method can capture vegetation extreme events and that comparison on a global scale is feasible. More importantly, agreement between the patterns generated by different methods covering different time spans indicates a consistent regularity and suggests that evaluation of vegetation vulnerability using vegetation extreme frequency is reasonable.

Relationships between vegetation extremes and climate variables
The results of this research have also shown that regionally differentiated patterns of vegetation extremes are closely associated with climate patterns, with vegetation in semi-arid and semi-humid regions being the most vulnerable (figure 3). Except for Cfc (oceanic) climates, which occur only at the southern end of Chile and in France and are adjacent to Cfb (maritime temperate climates, such as parts of Spain and France), climate types with high average vegetation extreme frequency are situated in sub-humid and sub-arid regions where vegetation is affected by drought, including Csa (dry-summer subtropical), Csb (Mediterranean), Dfb (warm-summer continental), BWh (low-latitude desert), BSk (middle-latitude steppe), Cwb (dry-winter temperate highland tropical), and Aw (savanna) (figure 3). The global pattern of vegetation extreme frequency also highlights the Amazon region, with the most severely affected areas in the eastern part, which has Am (tropical monsoon) and Aw (savanna) climates and could also be affected by water deficiency.
Unlike extremes in temperature and PDSI, extremes in precipitation showed a significant relationship with vegetation extremes in all biomes, indicating that precipitation is the major controlling variable of vegetation extremes on a global scale. Although PDSI takes the effects of both temperature and precipitation into account and has been very useful in detecting vegetation drought (Dai et al 2004), it showed a poor relationship with vegetation extremes in this study. This might be partially attributable to a mismatch in spatial resolution, because PDSI is the only variable with a resolution of 2.5 • ; however, this explanation needs further investigation and discussion. In the linear regressions between precipitation extreme frequency and vegetation extreme frequency, the slope increased with τ in quantile regression for most biomes, which agrees with the assumption that ignored microclimates could work as refuges and reduce the sensitivity of certain pixels. The highest slopes were found in temperate broadleaf forest and temperate grassland, suggesting high sensitivity to precipitation extremes. This result corresponds to the vegetation extreme frequency pattern that highlighted the vulnerability of semi-humid and semi-arid regions. If future precipitation extremes become more severe, these regions are likely to have more unstable vegetation or even to experience vegetation shifts.
The results of this research contradict the assertion that semi-arid and semi-humid regions, lacking the rapid response to drought displayed by other regions, might be more resistant to water deficits (Vicente-Serrano et al 2013). By finding the highest correlation between vegetation indices and the Standardized Precipitation Evapotranspiration Index (SPEI), Vicente-Serrano et al (2013) evaluated the drought-response time scales of different biomes. Their results generally agree with the present research on the slight impact of drought on cold and wet biomes (except for Amazonia in this study). However, it has been argued here that long vegetation-response time lags do not necessarily indicate resistance to drought, as suggested by Vicente-Serrano et al (2013), but are equally likely to result from more lasting or even more profound impacts. In addition, although vegetation in these regions has become physiologically adapted to regular periods of water deficit and might be resistant to non-severe water deficits, it might still display high sensitivity to severe water deficits, as in the precipitation extremes in this study. Many semi-arid and semi-humid regions are situated in ecotones and include the environmental limitations of different vegetation types, where drought or other environmentally adverse effects could have significant ecological consequences. The transition between savanna and forest, for example, is sensitive to drought (Staver et al 2011). Nevertheless, as Vicente-Serrano et al (2013) pointed out, determining drought vulnerability requires that the recovery time after drought be examined. The discrepancy between the results of these two studies might be explained by the difference in drought severity once recovery has been examined, which can be quantified only after the vegetation recovery time.
The global pattern of ecosystem vulnerability under trends of observed climate change overlaps extensively with the present results, especially for the Amazon rainforest, southeastern Africa, Central America, and the Indochinese peninsula, which indicates a strong relationship between climate change and vegetation decline (Gonzalez et al 2010). However, because the patterns acquired in this research were based on vegetation activity, these two patterns also exhibit inconsistencies which illustrate the complexity of this relationship. Certain regions that were highlighted only in this study require more attention because they were not captured in the climate changes described by Gonzalez et al (2010). Present predictions of vegetation change are largely based on observed climate change and climate projections together, and therefore these regions might display unexpected behavior in the future. They include northwestern and southeastern Australia, the central and eastern United States, and tropical and subtropical dry forest in southern Africa (Angola, Zambia, Zimbabwe, and Mozambique).

Uncertainties and prospects of this study
The major uncertainties in this study came mainly from the data sources. Although the GIMMS-NDVI data set has been corrected to improve its indication of vegetation condition, non-vegetation information still exists in it (Kobayashi and Dye 2005). The resampling technique used in this study could avoid interference only within each 0.5 • longitude × 0.5 • latitude grid cell, and therefore the results were still affected by persistent clouds and aerosols that could lead to data bias, particularly in tropical regions Dye 2005). The apparent leaf flush in the dry season and the LAI (leaf area index) in the wet season in the Amazon rainforest might not be explicitly associated with vegetation activity (Asner and Alencar 2010), which adds an extra element of uncertainty in detecting vegetation dynamics though NDVI. In addition, although NDVI has long been used successfully to monitor vegetation condition in low-vegetation-cover regions such as deserts (Dall'Olmo and Karnieli 2002, Richard and Poccard 1998, Tucker et al 1991, NDVI in these regions might be less accurate due to the possible influence of non-vegetation land cover (Kogan and Zhu 2001).
The relatively small sample size (25 years for each pixel for each two-week period) for extreme event detection, especially in the Box-Cox transformation, which lacks a method for quantifying the uncertainty of estimating λ, also increased random error. This uncertainty could be reduced by using a data set acquired over a longer period. Although NDVI values derived from MODIS data are available for 2006 to present, only the GIMMS-NDVI data set was used here because of possible inconsistencies introduced by a change in data source. For convenience in comparison, global vegetation was divided into eight biomes, with each biome made up of different vegetation types with different responses to climate change. Although linear models capture the overall relationships, they dampen the variation in characteristics within each biome. This is best shown in tropical and subtropical grassland and forest, where the diversity of vegetation types is high. Although the slopes of the linear models are significant, the discrete data distribution shows a mixture of vegetation with different sensitivities, which will require further investigation at a finer scale in the future.

Conclusions
Vegetation extreme frequency exhibits a pattern of global vegetation vulnerability that is consistent with earlier studies using different methods over different time spans. This pattern highlights highly vulnerable regions: the Sahel, tree savanna and temperate steppe in southern Africa, forest-prairie ecotones in North America, Central America, Amazonia, southeastern and northwestern Australia, and Mediterranean vegetation in southern Europe. Generally, vulnerability is highest in the Amazon rainforest, semi-arid regions, and semi-humid regions, specifically at MAT of 10-28 • C and MAP of 300-1000 mm, and is low in cold regions.
Compared to temperature and PDSI, extremes in monthly precipitation have a more significant effect on vegetation extremes at a global scale. The highest slopes in linear correlation analysis occur in temperate broadleaf forest and temperate grassland, suggesting that these environments might be the most sensitive to precipitation extremes among the eight biomes analyzed. This result confirms the high vulnerability of semi-arid and semi-humid regions and indicates that vegetation in semi-arid and semi-humid regions is more likely to be endangered under future climate change, especially changes in precipitation.