Remote Sensing of Environment 187 (2016) 14–29

Contents lists available at ScienceDirect

Remote Sensing of Environment journal homepage: www.elsevier.com/locate/rse

Diverse relationships between forest growth and the Normalized Difference Vegetation Index at a global scale Sergio M. Vicente-Serrano a,⁎, J. Julio Camarero a, José M. Olano b,c, Natalia Martín-Hernández a, Marina Peña-Gallardo a, Miquel Tomás-Burguera d, Antonio Gazol a, Cesar Azorin-Molina e, Upasana Bhuyan f, Ahmed El Kenawy g a

Instituto Pirenaico de Ecología, Consejo Superior de Investigaciones Científicas (IPE–CSIC), Zaragoza, Spain Área de Botánica, Departamento de Ciencias Agroforestales, EU de Ingenierías Agrarias, Universidad de Valladolid, Campus Duques de Soria, 42004 Soria, Spain Sustainable Forest Management Research Institute, Universidad de Valladolid and INIA, Avda. de Madrid 44, 34004 Palencia, Spain d Estación Experimental de Aula Dei (EEAD-CSIC), Zaragoza, Spain e Department of Earth Sciences - Regional Climate Group, University of Gothenburg, Sweden f Technische Universität München, Freising, Germany g Department of Geography, Mansoura University, Mansoura, Egypt b c

a r t i c l e

i n f o

Article history: Received 1 March 2016 Received in revised form 9 September 2016 Accepted 3 October 2016 Available online xxxx Keywords: Dendrochronology Normalized Difference Vegetation Index (NDVI) Tree rings Forest growth

a b s t r a c t This study compared the densest available database of tree-ring growth with the longest Normalized Difference Vegetation Index (NDVI) information available at the global scale to quantify the relationship between annual forest growth and the NDVI across different forest types and regions and to characterize the patterns of response of forest growth to NDVI values at different temporal scales. We found a general positive relationship between the inter-annual NDVI variability and the annual tree growth in most of the analyzed forests. Nevertheless, there were strong differences in the tree growth responses to NDVI, given that the annual tree-ring records in each forest responded in a different way to the magnitude, seasonality and accumulation period of the NDVI values. Thus, we found eight main patterns of tree-ring response to the NDVI, which were related to the forest type and climate conditions of each corresponding site. The identified patterns may be useful for determining early-warning signals of changes in forest growth over large areas based on remote sensing information. © 2016 Elsevier Inc. All rights reserved.

1. Introduction Forest productivity and growth are proxies of ecosystem health (Dobbertin, 2005). They show tight relations with climate variability (Barber et al., 2008; Lloyd and Fastie, 2002; Hirota et al., 2011), including extreme events, such as droughts (Martínez-Vilalta and Piñol, 2002; Bigler et al., 2006; Pasho et al., 2011). In a warmer climate, forest dieback and related tree mortality episodes are often associated with intensified drought stress (Phdersen, 1998; Williams et al., 2010, 2013; Allen et al., 2010; Liu et al., 2013). In this regard, assessing forest growth is usually carried out by national administrations or research groups, which develop permanent plots based inventories in order to record changes in stem diameter over large areas (Fang et al., 2001; Fridman and Walheim, 2000; Jenkins et al., 2003). While information on forest growth using this stand scale is highly detailed, it is extremely expensive due to the intensive field sampling. Also, the temporal resolution

⁎ Corresponding author. E-mail address: [email protected] (S.M. Vicente-Serrano).

http://dx.doi.org/10.1016/j.rse.2016.10.001 0034-4257/© 2016 Elsevier Inc. All rights reserved.

of this information is very low since these inventories are usually repeated every five years. The strong inter-annual variability that characterizes the current climate scenarios (Hartmann et al., 2013) makes it necessary to assess forest productivity and growth to provide a rapid feedback for management and conservation practices. The most common approach to measure year-to-year variability in radial growth is to use tree-ring records, such as width (Fritts, 2001). Tree-ring samples allow for determining how much wood is formed by the stem, which can be directly related to biomass gain and carbon uptake (Ketterings et al., 2001; Van Breugel et al., 2011; Babst et al., 2014a, 2014b). This approach is widely used in forest management and ecological research, including a profuse research body that focuses on how climate change may affect forest growth in the recent decades (Peterson and Peterson, 2001; Tardif et al., 2003; Sarris et al., 2007; Williams et al., 2010, 2013, Shestakova et al., 2016). However, although tree rings provide an accurate and retrospective measure of productivity and growth, their use also involves intensive field work to collect wood samples, in addition to time-consuming laboratory work to cross-date and measure the collected material (Cook and Kairiukstis, 1990). These drawbacks limit the

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

potential to use this approach not just for monitoring real-time forest growth over large spatial scales, but also for management objectives, such as getting early-warning signals of forest dieback (Camarero et al., 2015). Remote sensing images provide spectral radiance in different regions of the electromagnetic spectrum, offering the possibility to obtain indices that inform of the vegetation conditions based on the differential absorption by vegetation (Rouse et al., 1974; Tucker, 1979). In the past four decades, information provided from earth observation satellites has gained wide use for several ecological and forestry applications, including- among others: (i) the quantification of the impact of drought events on forest productivity (Vicente-Serrano, 2007; Büntgen et al., 2010); (ii) the characterization of succession processes (Hall et al., 1991; Song and Woodcook, 2003); (iii) the determination of forest biomass (Lu, 2006; Timothy et al., 2016); (iv) the assessment of tree density (Crowther et al., 2015); and (v) the management of forest resources (Takao et al., 2010). Thus, intense efforts have been devoted to determining forest biophysical variables by means of remote sensing data. The most popular vegetation index is the Normalized Difference Vegetation Index (NDVI), which measures the fractional absorbed photosynthetically active radiation (Myneni et al., 1995) and exhibits a strong relationship with vegetation parameters, such as Leaf Area Index (LAI) (Baret and Guyot, 1991, Carlson and Ripley, 1997; Shabanov et al., 2005), green biomass (Gutman, 1991, Cihlar et al., 1991, Wylie et al., 2002), and vegetation cover (Duncan et al., 1993; Gillies et al., 1997; Gouveia et al., 2016). The NDVI properties have allowed the use of this information for estimating the Net Primary Production (NPP) (Goward and Dye, 1987; Running et al., 2004; Hasenauer et al., 2012). Different studies have already found a strong relationship between NPP and radial growth (e.g., Granier et al., 2008; Babst et al., 2013, 2014a, 2014b; Vicente-Serrano et al., 2015), albeit with significant differences, particularly those related to species, sites and environmental conditions (Rossi et al., 2006). While the available information from different satellite platforms could provide an opportunity to estimate the tree growth over large forested areas, there are very few approaches that have tested the use of either primary remote sensing information (i.e., reflectance values) or vegetation indices, as related to the inter-annual variability of forest growth (e.g., Vaganov et al., 1999; Kaufmann et al., 2004, 2008; Beck et al., 2011; Berner et al., 2013). This shortcoming is mainly associated with the lack of long time series of satellite imagery. The continuous remotely sensed data with adequate temporal resolution are available from the 1980s (NOAA-AVHRR satellites); other earth observation systems cover either shorter periods (b 15 years) (e.g., MODIS or SPOT-vegetation) or being affected by large temporal data gaps (e.g., Landsat). Importantly, the available spatial resolution of the satellite data is too coarse to be related directly to growth data from small forest areas, which makes the availability of suitable treering information an important constraint. Albeit the strong economic and environmental implications of estimating the temporal variability of forest growth at large spatial scale, studies that relate forest growth (e.g., tree-ring width) and satellite vegetation indices at a global scale are still lacking. Here we compare the densest available database of tree-ring growth with the longest NDVI information available at the global scale. Thus, the aims of this study are: (i) to quantify the relationship between annual forest growth and the NDVI across different forest types and regions worldwide; and (ii) to characterize the patterns of response of forest growth to NDVI values at different temporal scales, under different forest types and/or environmental conditions. To our knowledge, this is the first comprehensive study about this issue, covering different tree species and phenologies, from a global perspective.

15

2. Datasets 2.1. NDVI data Some global datasets, such as SPOT (Satellite Pour l'Observation de la Terre) Vegetation and MODIS (Moderate-Resolution Imaging Spectroradiometer) – Aqua, provide NDVI at high spatial resolution (1 km). Nonetheless, the temporal coverage of these data is too short (starting mostly in the early 2000s) for comparisons with the available longer records of tree-ring data. Alternatively, the unique earth observation program, which enables the analysis of vegetation activity over the past three decades, is provided by the National Oceanic and Atmospheric Administration (NOAA) Polar orbiting satellites. This platform employs the Advanced Very High Resolution Radiometer (AVHRR) sensor to collect the necessary spectral information for calculating NDVI. The most widely used global dataset of NDVI from NOAA-AVHRR data is the Global Inventory Monitoring and Mapping Studies (GIMMS) database (Tucker et al., 2005), which covers the period 1981–2013. The quality and consistency of the GIMMS data have already been assured by corrections for (i) sensor degradation, (ii) sensor inter-calibration differences, (iii) solar zenith and viewing angles, (iv) volcanic aerosols, (v) atmospheric water vapor, and (vi) cloud cover (Pinzon and Tucker, 2014). In this study, we used the GIMMS NDVI3g (http://ecocast.arc. nasa.gov/data/pub/gimms/3g.v0/; last access on 1st March 2016) provided at a spatial resolution of 0.083° (which corresponds to 9.23 km at the equator, 7.06 km at 40°N and 4.06 km at 60°N), with a biweekly temporal interval for the period from July 1981 to December 2013 (Pinzon and Tucker, 2014). For validation purposes, we have used the Monthly L3 Global (MOD13A3) NDVI dataset from the MODIS satellites (https://lpdaac. usgs.gov/dataset_discovery/modis/modis_products_table/mod13a3). For this purpose, we have selected the scenes corresponding to the location of each one of the available tree-ring series (see Section 2.2), and used the period between February 2000 and December 2013. This dataset covers a shorter period than the GIMMS dataset but it has higher spatial resolution (1 km), which allows to determine if the patterns of tree-ring/NDVI correlations obtained at low spatial resolution (GIMMS) are comparable to those obtained at more detailed spatial resolutions (1 km from the MODIS dataset). 2.2. Tree-ring data Here we used the tree-ring data available at the International TreeRing Data Bank (ITRDB; http://web.utk.edu/~grissino/itrdb.htm; last access on 1st March 2016) (Grissino-Mayer and Fritts, 1997). This repository contains tree-ring data collected from N2000 sites, covering six continents and representing N 100 tree species. These series are available online via http://www.ncdc.noaa.gov/paleo/treering.html; last access on 1st March 2016). Here, we restricted our analysis to tree-ring width series, which have at least 18 years of common data with NDVI data over the period 1981–2013. Following this criterion, a total of 702 site chronologies were defined globally. Fig. 1 shows the location of the tree-ring sites with time series available for this study. They mostly cover regions of North America, Europe and Asia. Each local chronology represents the average growth series of several trees (at least ten) of the same species growing at the same site. The wood samples were taken and processed following a standard protocol and taking, in general, two radial cores per tree at 1.3 m following Stokes and Smiley (1968). Tree-ring width is a good proxy of tree lifetime growth and carbon uptake as woody pools (Fritts, 2001). However, given that the ring width declines along a stem radius, as a function of time (age), it is necessary to convert raw tree-ring width into indexed (detrended) data (Cook and Kairiukstis, 1990). To accomplish this task, tree-ring width measurements provided to the ITRDB were standardized and detrended using ARSTAN software (Cook and Krusic, 2005). Individual series of tree-ring widths were fitted with negative

16

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Fig. 1. Spatial distribution of the 702 sites with available tree-ring series. Black points: 497 selected forests for analysis; white points: 205 non-selected forests. The color scale represents the mean NDVI across the world. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

exponential or linear functions and residuals were obtained by dividing the observed values by the fitted ones. For each site, the resulting residuals were then subjected to an autoregressive modelling and averaged for each year using a biweight robust mean to obtain a mean chronology of ring-width indices. Finally, the low- to medium-frequency variability and the first-order autocorrelation were removed to obtain the residuals or pre-whitened ring-width indices. 3. Methods 3.1. Selection of the spatio-temporal homogeneous forests To account for the relatively coarse spatial resolution of the satellite imagery used in this study, it was necessary to make an accurate selection of the available sites. The aim was to assure that the NDVI data are representative of the tree-ring growth at the forest scale. Therefore, we selected those sites that represent homogenous surrounding areas, with low spatial variability of NDVI values. Specifically, we selected the NDVI pixels at which tree-ring series were available. Then, we calculated the average values of NDVI in these pixels, and calculated the coefficients

of variation with surrounding pixels, considering a spatial window of 3 × 3 pixels. Fig. 2 shows the probability density function of the obtained coefficients of variation corresponding to each one of the 702 sites. In the majority of cases (94%), differences in the coefficient of variation between each particular NDVI pixel and the surrounding pixels were lower than 0.4. This reveals a similarity between the average NDVI in the forest sites and those of surrounding areas, suggesting a low landscape variability around each forest site. Considering a threshold of 0.4 as the difference in the coefficient of variation between each selected pixel and the surrounding pixels, we defined and removed the forests with high NDVI spatial heterogeneity. We also considered an additional criterion to retain those forests with the homogeneous temporal variability of NDVI in relation to the surrounding pixels. Specifically, we calculated the Pearson's r correlation coefficient between the biweekly NDVI series in the forest site and the NDVI series in a 3 × 3 pixel window. The aim was to determine whether the forest site showed a similar NDVI evolution, with respect to those of the surrounding pixels. The probability density function of the correlations metrices is summarized in Fig. 3. As illustrated, in 87% of the forest pixels, the correlation with nearby pixels was generally higher

Fig. 2. Relative frequency of the spatial coefficients of variation applied to the average NDVI values considering a spatial filtering of 3 × 3 pixels.

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

17

Fig. 3. Relative frequency of the Pearson's r correlation coefficients between the NDVI series corresponding to the forest site and those of the surrounding sites using a 3 × 3 window.

than 0.7, suggesting high temporal consistency of NDVI between the forest sites and the surrounding landscape. To reduce uncertainties associated with the temporal variability of the NDVI series, we selected only the forests whose temporal variability of NDVI showed high agreement with those of the surrounding pixels using a threshold of r = 0.7 (all the surrounding pixels may have a r value over this threshold to select the forest series). Following our procedure and the defined criteria, only 497 forests were chosen from a bulk of 702 forests. With this selection, we removed the forests located in highly fragmented landscapes or characterized by a very different NDVI temporal variability, as compared with the surrounding areas. This approach guaranteed that the selected forests were characterized by homogeneous NDVI values from a spatiotemporal perspective, which reduces the possible impacts of the coarse spatial resolution of NDVI on the final results. To test whether our criteria were sufficient to identify relatively homogeneous forests that represent the temporal variability of the NDVI values of the pixels in which they are located, we used a high resolution land cover map obtained from the MERIS sensor (Glocover, http://due. esrin.esa.int/page_globcover.php). This land cover map shows a good classification accuracy at the global scale (Bontemps et al., 2013). We calculated the surface percentage of the forest categories at each pixel that can ascribe the tree species from which tree-ring data were sampled. We found that the most representative categories were: i) closed needle leaved evergreen forests (15.4%), ii) closed broad leaved deciduous forests (11.5%) and iii) closed to open mixed broad leaved deciduous and evergreen forests (5.4%). The total percentage of the categories that correspond to different forest types was 80.7% in the selected 497 forests. Among them, the majority of the analyzed NDVI pixels (398) showed a percentage higher than 80% of the forest coverage, while only 47 cases exhibited a percentage lower than 40%. 3.2. Relationship between the NDVI and the tree-ring series We calculated Pearson's correlation coefficients between the annual tree-ring width indices and the de-trended biweekly NDVI series at each forest site. Since the cumulative Net Primary Production (NPP) over long periods can be a better estimator of tree-ring width than NPP of shorter periods (Gough et al., 2008; Zweifel et al., 2010), we correlated the annual tree-ring data with the NDVI summarized at different time-scales. For NDVI time scales, we referred to the average NDVI over the previous n biweekly periods. Therefore, we correlated the tree-ring data series with the 24 NDVI biweekly (i.e., two per month) series at

time-scales varying from 1 to 48 biweekly periods (i.e., two years). We considered the NDVI values not only for the corresponding year, but for the previous year as well. This is mainly because tree-ring growth may be impacted by tree activity during the previous year (e.g., photosynthesis, carbohydrate synthesis and storage and bud development) (Fritts, 2001). For each tree-ring width series, we obtained 576 correlations (24 semi-monthly periods × 24 time-scales). This procedure allowed for determining whether the ring-width indices are linked more to the NDVI values of the previous and/or the corresponding year. Results are summarized by means of surface plots, which also indicate the threshold for significant correlations (p b 0.05). 3.3. Assessment of the robustness of the experimental design using MODIS NDVI at 1 km resolution We have used the MODIS NDVI data set to determine if the patterns of tree-ring/NDVI relationships obtained at the spatial resolution of the GIMMS data are coherent at higher spatial resolutions (1 km). This assessment has required the degradation of the spatial resolution of the MODIS dataset at 2 × 2 km, 4 × 4 km and 8 × 8 km in the location of each 702 tree-ring samples by averaging the NDVI values. Then we analyzed the relationship between the 1 km NDVI at the location of the tree-ring samples and the 2-, 4- and 8-km NDVI data. In this analysis we computed the correlation coefficient considering (i) the entire time-series between 2000 and 2013 in form of deseasonalized anomalies (Fig. 4A) and (ii) the 12 monthly time-series (Figs. 4B, C and D). As expected, the correlation between the 1 km NDVI series and the NDVI series at lower spatial resolutions shows a decrease in the agreement as the spatial resolution decreases. The series degraded at 2 × 2 km show high agreement with the original series at 1 km, whereas the correlation decreases at the spatial resolutions of 4- and 8-km. Nevertheless, the correlations are dominantly positive and statistically significant in the majority of the analyzed forests, even at the lowest spatial resolution considered (8 × 8 km). Thus, the average of correlation between the forest NDVI series at 1 km and at 8 km is higher than r = 0.6, independently of the month of the year. This means that low-resolution (4- to 8-km) NDVI series may be representative of the general interannual variability of the vegetation activity in the available forest sites. We have also selected the annual tree-ring data series with at least 8 years of common period with the MODIS NDVI dataset and correlated with the 1- to 24-month cumulative NDVI series. 155 series were available, 102 of them corresponding to the selected forests according to the

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

1.0

1.0

0.5

0.5

Pearson's r

Pearson's r

18

0.0

-0.5

0.0

-0.5

-1.0 2x2

4x4

8x8

1

1.0

1.0

0.5

0.5

Pearson's r

Pearson's r

2x2

B)

A) -1.0

0.0

2

3

4

5

6

7

8

9 10 11 12

0.0

-0.5

-0.5

4x4

C)

8x8

D) -1.0

-1.0 1

2

3

4

5

6

7

8

9 10 11 12

1

2

3

4

5

6

7

8

9 10 11 12

Fig. 4. Box-plots with the correlation coefficients between the time series of NDVI in the 702 available tree-ring series across the world at the spatial resolutions of 1 × 1-, 2 × 2-, 4 × 4- and 8 × 8-km. A) Correlations from the complete series of NDVI transformed to deseasonalized standard anomalies to allow comparisons. B, C and D) monthly series.

criteria detailed in Section 3.1. Fig. 5 shows an example of the correlation patterns in three of the available forests between the annual treering series and the monthly NDVI at different time-scales. In these three forests, although the magnitude of correlations may differ, the patterns of tree-ring/NDVI correlations for the different months and NDVI time-scales are similar, suggesting that the response of tree-ring is mostly insensitive to the spatial resolution of the NDVI series. Thus, this is the dominant pattern when analyzing all 155 forests. We have analyzed the agreement between the patterns of tree-ring/ NDVI correlations obtained at the spatial resolutions of 1- and 8-km in each one of the 155 forests. For this purpose we calculated the correlation between the patterns of tree-ring/NDVI relationships obtained from the MODIS NDVI series at 1-km and 8-km of spatial resolution. Fig. 6 shows relative histograms for the 155 forests but also for the selected (102) and non-selected (53) forests according to the criteria explained in Section 3.1. In the majority of the 155 forests analyzed, the pattern of the tree-ring/NDVI correlation obtained at 1-km is quite similar to that obtained with the MODIS NDVI series at 8-km. Therefore, the 47.8% of the 155 series shows correlations between both patterns higher than r = 0.75, and the 68% of the forests a correlation above r = 0.5. The relative frequencies are similar for the non-selected and selected forests, although the frequency of correlation values is higher in the selected than in non-selected forests (50.5% vs. 43.4% with correlations higher than 0.7, and 70.3% and 60% with correlations higher than 0.5, respectively). Finally, we have analyzed the magnitude of the maximum correlation between 1- to 24-month NDVI time-scales and the annual treering growth series with the NDVI series at 1- and 8-km. The maximum correlation obtained at both spatial resolutions shows high agreement in all the 155 forests (Fig. 7), being the agreement higher in the selected than in the non-selected forests. This pattern is also identified for the month in which the maximum correlation is found (r = 0.3, r = 0.19

and r = 0.35 for the all available, the non-selected and the selected forests, respectively) and also for the NDVI time-scale at which maximum correlation is found (r = 0.25, r = 0.18 and r = 0.29, respectively). All these results based on the high resolution MODIS-NDVI data reinforce the use of the relatively low-resolution GIMMS-NDVI data since the patterns of tree-ring/NDVI correlations obtained using low spatial resolution NDVI (8 km) are comparable with those obtained at 1 km of spatial resolution in the majority of the analyzed forests (Fig. 7). Although the robustness of the relationship has been found for both selected and non-selected forests according to the criteria explained in Section 3.1, we have focused our analysis on the 497 forests, in which the NDVI signal is representative of larger surrounding areas to avoid residual geocoding problems in the data. 3.4. Summarizing correlation patterns Variability in the correlations between the tree-ring growth series and the NDVI series was summarized using the Principal Component Analysis (PCA; Richman, 1986). Then, chronologies were classified on the basis of the similarities of correlations calculated between the tree-ring and the NDVI series at the different time-scales. In order to identify the general spatial patterns in the calculated correlations, this classification was made using an S-mode of PCA. We used a correlation matrix to calculate the PCA (Barry and Carleton, 2001). We obtained the components in the original correlation coefficient values using the weight coefficients of each forest in each component. The number of retained components was defined based on the percentage of the total explained variance, as suggested by the scree-plot. Given that the maximum correlation with any retained component could be positive or negative, we assigned the forests to different categories, based on the component in which the forest showed the strongest correlation coefficient, but taking into account the sign of the correlation (i.e., negative or

22

22

20

20

20

18

18

18

16

16

16

14 12 10

Time-scale

22

Time-scale

24

14 12 10

14 12 10

8

8

8

6

6

6

4

4

1x1

2 1

2

3

4

5

6

7

8

9

10

11

12

4

1x1

2 1

2

3

4

5

Month

6

7

8

9

10

11

1

12

22

22

20

20

20

18

18

18

16

16

16

14 12 10

Time-scale

22

14 12 10

8

6

6

8x8

4

5

6

7

Month

8

9

10

11

12

4

8x8

2 4

1

2

3

4

5

6

7

Month

6

7

8

9

10

11

12

8

9

10

11

8

9

10

11

12

8x8

2 12

-0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8

10

6

3

5

12

8

2

4

14

8

1

3

Month 24

2

2

Month 24

4

1x1

2

24

Time-scale

Time-scale

cana426 -118.03ºW 52.55ºN

germ160 8.59ºW 51.09ºN 24

1

2

3

4

5

6

7

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Time-scale

chin071 99.55ºW 38.12ºN 24

Month

Fig. 5. Pearson's r correlation coefficients calculated between the time series of tree-ring width and cumulative MODIS NDVI at different (1–24 months) temporal scales and two different spatial resolutions (1- and 8-km) for selected three forests. Dotted lines denote values, which are statistically significant at p b 0.05.

19

20

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Non-selected

All available

Selected

Relative frequency (%)

12 10 8 6 4 2 0 -1.0

-0.5

0.0

0.5

1.0

-1.0

-0.5

0.0

0.5

1.0

-1.0

-0.5

Pearson's r

Pearson's r

0.0

0.5

1.0

Pearson's r

Fig. 6. Relative frequencies of Pearson's r values between the patterns of MODIS-NDVI/tree-ring relationships at the NDVI spatial resolutions of 1- and 8-km. The frequency histograms for all available forests (155), and non-selected (53) and selected forests (102) are shown.

positive). We also mapped the loadings, which summarize the correlation between the tree-ring growth and NDVI for each particular component, as well as the general pattern (PC) that represents a number of forests. 3.5. Factors explaining the spatial differences in the NDVI-tree-ring correlations

Predictive Discriminant Analysis (PDA), which explains the value of a dependent categorical variable based on its relationship to one or more predictors (Huberty, 1994; Hair et al., 1998). PDA allowed for assessing which predictors contributed more to the intercategory differences of the PCs that summarized the tree-ring growth/NDVI dependency. 4. Results

We used several sources of information to assess the influence of different biophysical and climate conditions on the different responses of the tree-ring width chronologies to NDVI variations computed at different time-scales. First, we assessed the influence of leaf habit of each forest species (i.e., deciduous vs. evergreen species). Second, we focused on the role played by a range of climatic variables (e.g., annual precipitation and mean air temperature). In this regard, we also assessed the impact of water balance, defined as precipitation minus reference evapotranspiration, providing evidence on aridity intensity. Climate data were obtained from the CRU TS 3.21 data set (https://crudata.uea.ac.uk/cru/ data/hrg/; last access on 1st March 2016; Harris et al., 2014). Third, the possible influence of some common geographical variables (e.g., latitude and elevation) was also assessed. In order to summarize the responses of ring-width indices to the NDVI time-scales, average values for these geographical variables were obtained for the selected PC groups. The contribution of these explanatory variables to the spatial differences in the ring-width responses to NDVI at different timescales was estimated using a

All Available

4.1. Relationships between the tree-ring growth and NDVI Our findings revealed that almost 67% of the forests showed a positive and significant correlation (p b 0.05) with NDVI, regardless of the cumulative period or the month of the maximum correlation. However, there was a strong diversity in the patterns of the correlations between the tree-ring width indices and the NDVI calculated at different timescales. Fig. 8 shows some examples of representative forest types and locations. Although there were very different patterns of response of forest growth to the NDVI, growth is generally more related to the cumulative NDVI values over long periods than to the NDVI recorded over a particular biweekly period. For example, in a Spruce (Picea spp.) forest located in Damra (Bhutan), the annual tree-ring growth is mainly explained by the cumulative NDVI during three biweekly periods in June. In a Blue oak (Quercus douglasii) forest located in the American Canyon (California, USA), forest growth was principally linked to the cumulative NDVI values from March to May, although

Non-selected

Selected

Maximum correlation (8 x 8)

1.0

Pearson's r = 0.59

Pearson's r = 0.53

Pearson's r = 0.61

0.8

0.6

0.4

0.2

0.0 0.0

0.2

0.4

0.6

0.8

Maximum correlation (1 x 1)

1.0

0.2

0.4

0.6

0.8

Maximum correlation (1 x 1)

1.0

0.0

0.2

0.4

0.6

0.8

1.0

Maximum correlation (1 x 1)

Fig. 7. Relationship between the magnitude of the maximum correlation between 1- to 24-month NDVI time-scales and the annual tree-ring growth using the MODIS-NDVI series at 1- and 8-km resolution. The analysis is shown for the all available forests but also for the non-selected and selected forests.

Spruce, Damra (Buthan), 27.28ºN-89.31ºE

45 40 35 30 25 20 15 10 5 1

2

3

4

5

6

7

8

9

45 40 35 30 25 20 15 10 5 1

10 11 12

2

3

4

5

40 35 30 25 20 15 10 5 3

4

5

6

7

8

9

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

Qilianshan Juniper, Delingha (China), 37.27ºN-97.32ºE

45

2

35 30 25 20 15 10 5 5

6

7

10 11 12

40

-0.6

35

-0.4

30 -0.2

25 20

0.0

15

0.2

10 0.4

5

0.6

1

8

9

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

40

4

9

2

3

4

5

6

7

8

9

10 11 12

Month

Scots Pine, Nordrhein-Westfalen (Germany), 50.35ºN-6.28ºE

3

8

Lodgepole Pine, Fourth of July Hill (USA), 37.52ºN- 119.22ºW

10 11 12

45

2

7

45

Month

1

6

Month

Month

1

21

Blue Oak, American canyon (USA), 35.18ºN- 120.16ºW

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Deodar Cedar, Zairat Chitral (Pakistan), 35.21ºN-71.48ºE

45 40 35 30 25 20 15 10 5

10 11 12

1

2

3

4

Month

5

6

7

8

9

10 11 12

Month

Fig. 8. Pearson's r correlation coefficients calculated between the time series of tree-ring width and cumulative NDVI at different (1–48) temporal scales for selected six forests in the world. Black points correspond to the maximum correlation found for each forest. Dotted lines denote values, which are statistically significant at p b 0.05.

100 90 80 70

% of variance

longer NDVI cumulative periods in the following months (July– September) were also strongly linked to the tree-ring growth. These patterns were very different from those observed in the Lodgepole pine (Pinus contorta) forest situated in the Fourth of July Hill (California, USA), where the maximum correlations between tree growth and NDVI were recorded at the cumulative 33-biweekly NDVI values prior to May. This means that the recorded NDVI over the main part of the previous year had a positive influence on tree growth in the following year. We also found a wide range of patterns, as observed in other forests, such as Delingha (Central China), the temperate Scots pine (Pinus sylvestris) forests (Germany) and Zairat Chitral (NW-Pakistan). We also obtained the main patterns of growth–NDVI correlations from the 497 selected forests by applying a Principal Component Analysis (PCA). Scree-plot criterion suggested keeping the first four PCA axes, which summarized the tree-ring growth/NDVI dependency and explained 70% of the total variance (Fig. 9). These components reflected the variability in the patterns of response of tree growth to cumulative NDVI values. Fig. 10 illustrates the spatial distribution of the PC loadings corresponding to the selected four components, while Fig. 11 shows the patterns of the four retained PCs, considering both positive and negative correlations. Component 1 (+) is characterized by positive and

60 50 40 30 20 10 0 0

5

10

15

20

25

30

Components Fig. 9. Total explained variance (%: bars) and the cumulative variance (%: black line) obtained following the Principal Component Analysis.

22

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Fig. 10. Spatial distribution of the loadings of the four principal components, summarizing the patterns of correlation between the tree-ring width indices and NDVI at a global scale.

significant correlations between growth and cumulative NDVI for 3–4 biweekly periods between July and August. It can be noted that growth did not respond to the NDVI values during other months of the year and cumulative periods. While the location of the forests represented by this pattern is not clear, given the high spatial variability in the magnitude of the correlation coefficients at the global scale (Figs. 10 and 11), this pattern can generally be seen in the forests of the Mediterranean region (mainly the Iberian Peninsula and Turkey), central and western areas of North America and central Asia. In contrast, component 1 (−) revealed a very different pattern of response since the tree-ring growth showed better correlation with the cumulative NDVI in February, considering the previous 35 biweekly periods as the best accumulation period of the NDVI. This finding suggests that the NDVI over the previous year is more important than NDVI of the current year in explaining the inter-annual variability of tree growth in forests corresponding to this pattern (e.g., the forests of Eastern North-America and Tasmania Islands). Component 2 (+) indicated high correlations with the cumulative NDVI recorded at long-time scales, with a maximum correlation recorded at the 45 biweekly periods in July, suggesting that the accumulation of NDVI over the previous and current year until September showed a stronger predictive capacity of tree growth. This pattern was mainly found in the forests of central and South Asia, South Canada, North USA and Siberia. On the contrary, component 2 (−) was characterized by a poor response of the tree growth to NDVI variability, which was only recorded in particular cold regions with boreal summers. Specifically, this pattern corresponded to the forests from cold biomes, such as those located in the Alps mountains or in the boreal forests in Northern Canada and Scandinavia. In these regions, NDVI values are typically null over large periods of the year, as a consequence of low air temperatures and snow coverage. Component 3 (+) suggested a clear pattern of growth response to the cumulative NDVI values between February and June, with a maximum correlation found at the 5 biweekly periods in May. Although there were no clear spatial patterns of this component, it was recorded in some forests of East USA, central Europe, North Mongolia and Chile. For component 3 (−), a negative correlation dominated during the first half of the year considering the cumulative periods between 1 and 7 months, though being statistically

insignificant (p N 0.05). In contrast, significant positive correlations were recorded in the late summer and early autumn. This pattern was represented in the forests of North India and Pakistan, besides the forests of North America. Finally, PC4 exhibited a clear opposite pattern between the positive and negative patterns. While the positive pattern showed negative correlations considering a range of cumulative NDVI periods and months of the year, the negative pattern revealed positive correlations between tree growth and the cumulative NDVI over long periods. These two contrasted patterns appeared in a small number of forests in different regions worldwide. Fig. 12 summarizes some statistics corresponding to the forests classified in each one of the eight groups: the magnitude of maximum correlation, the month (considering biweekly periods) in which maximum correlation was recorded and the timescale at which the NDVI values were accumulated to reach the maximum correlation. Component 1 (+) and component 1 (−) showed similar maximum correlations between NDVI and tree growth, albeit with large differences in terms of the month and the timescale at which the maximum correlation was reached. For component 1 (+), the maximum correlation was noted in July, while the maximum correlation corresponding to the negative mode was observed in March considering long time scales. For component 1, there were no statistically significant differences in the magnitude of correlations, whereas significant differences were recorded considering the biweekly period and the timescale of the correlation. As noted, the maximum correlation among the different groups was found for component 2 (+). This correlation was dominantly recorded at the beginning of the summer (similar to component 1 (+), but at longer timescales). In contrast, component 2 (−) showed low correlation values. As observed for component 1, the two modes of component 3 revealed similar correlation values, but recorded in different periods: mid May for the positive and September for the negative. Both modes considered time scales between 2.5 and 3 months. Finally, component 4 (+) showed low correlation values, which corresponded to short timescales. In contrast, the negative mode was characterized by high correlation coefficients, dominantly recorded in March and considered the cumulative NDVI values over the previous three months.

Component 1 (positive) 45 40 35 30 25 20 15 10 5 1

2

3

4

5

6

7

8

9

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Component 1 (negative) 45 40 35 30 25 20 15 10 5 1

10 11 12

2

3

4

5

45 40 35 30 25 20 15 10 5 3

4

5

6

7

8

9

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

Component 2 (positive)

2

40 35 30 25 20 15 10 5 4

5

6

7

8

9

-0.7

30 -0.6

25 20

-0.5

15

-0.4

10

-0.3

5 -0.2

2

3

4

5

35 30 25 20 15 10 5 5

6

7

8

9

10 11 12

-0.1

Month

0.0

Component (negative) Component 33 (negative)

0.1 0.2 0.3

35 35 30 30

0.4

25 25 20 20

0.6

0.5

0.7

15 15 10 10 55 1

8

9

Cummulative Bi-weekly NDVI periods

Cummulative Bi-weekly NDVI periods

40

4

7

2

3

4

55

66

77

88

99 10 10 11 11 12 12

Month Month

Component 4 (positive)

3

6

45 45 40 40

10 11 12

45

2

10 11 12

35

Month

1

9

40

1

Cummulative Bi-weekly NDVI NDVIperiods periods Cummulative Bi-weekly

Cummulative Bi-weekly NDVI periods

45

3

8

45

10 11 12

Component 3 (positive)

2

7

Component 2 (negative)

Month

1

6

Month

Month

1

23

10 11 12

Month

Component 4 (negative) 45 40 35 30 25 20 15 10 5 1

2

3

4

5

6

7

8

9

10 11 12

Month

Fig. 11. Main modes of correlation between tree-ring width indices and 1- to 48-long cumulated biweekly NDVI values (y-axis).

4.2. Factors explaining the different growth responses to NDVI Fig. 13 depicts the average values and standard errors of the annual air temperature, precipitation, climatic water balance, latitude and elevation corresponding to the eight obtained patterns of NDVI-growth correlations. The forests characterized by an NDVI-growth correlation that resembles component 2 (−) and component 4 (−) are generally

located in the coldest regions, while forests corresponding to component 4 (−) are situated in more humid areas and at lower latitudes. Component 1 (+) and component 1(−) did not show significant differences in the annual mean air temperature, but there were contrasting differences in the annual mean precipitation and water balance, since forests associated with component 1 (+) are generally located in drier sites than those of component 1 (−). There were no statistically

24

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Fig. 12. (Left) Mean values of the maximum correlation between the tree-ring growth and the NDVI; (central) the biweekly period at which the maximum correlation is recorded; and (right) the time-scale (the cumulative NDVI period) that correspond to the maximum correlation for the correlation patterns between tree-ring growth and NDVI obtained by means of the Principal Component Analysis. Error bars indicate two standard errors of the average.

significant differences in the latitude and elevation at which the forests of component 1 (+) and component 1 (−) are located. Forests of component 2 (+) are mainly assigned to areas with low precipitation (b600 mm/year). In these regions, the magnitude of the atmospheric evaporative demand is low, leading to the least negative water balance among all forests. Differences in the location of component 3 (+) and component 3 (−) patterns are largely determined by annual precipitation variations. As illustrated, component 3 (−) forests are mainly located in areas that record lower precipitation than those of component 3 (+). Finally, we noted that forests of component 4 (+) are located at the lowest altitudes, while forests corresponding to component 4 (−) are associated with the most lowest latitudes, relative to all other components. The Predictive Discriminant Analysis (PDA) allowed for determining the relative importance of the various climatic and environmental factors in assessing the relationship between NDVI and tree growth. Our results revealed three significant functions, which accounted for 70.9% of the total variance. The first function accounted for 36.5% of the total variance, with climatic water balance and precipitation being the variables with the most important correlation with the discriminant function (Table 1). This function mainly explains the water availability during the boreal summer and also at the annual scale. The coefficients corresponding to the variables of this function showed positive values over the different seasons, with the exception of precipitation and

water balance during the boreal winter. Evergreen and deciduous forests showed contrasted patterns, with positive and negative coefficients, respectively. The second function accounted for 20.3% of the total variance. The contrast between evergreen and deciduous forests was noted again and a negative coefficient corresponding to the annual water balance was observed. Finally, the third function represents precipitation and water balance during the boreal spring and autumn and annually, accounting for 14.1% of the total variance. Function 1 clearly separated between component 1 (+) and component 1 (−), but also between component 4 (+) and component 4 (−) (Fig. 14). Component 1 (+) and component 4 (+) showed negative values for function 1. Forests characterized by these growth responses to NDVI correspond mainly to deciduous forests located in areas affected by water deficit during the boreal summer. The opposite pattern characterizes component 1 (−) and component 4 (−), which corresponds to evergreen forests located in more humid areas. The remaining groups related to components 2 and 3 are assigned to evergreen forests located in areas characterized by high humidity during the boreal summer. Nevertheless, discriminant function 1 did not establish clear differences between component 2 and component 3 groups. Component 1 (+) and component 4 (+) were differentiated by function 2. This reveals that forests of component 4 (+) are located in drier areas than those of component 1 (+). Also, higher differences between component 2 and component 3 groups emerged from the analysis of discriminant

Fig. 13. Mean values of the annual mean air temperature, precipitation, climatic water balance, latitude and elevation of the different correlation patterns obtained by means of the Principal Component Analysis. Error bars indicate two standard errors of the average.

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29 Table 1 Structure matrix of the first three discriminant functions obtained from the Predictive Discriminant Analysis (PDA). The table shows the correlation values of each predictor variable with the three discriminant functions. Abbreviations: DJF, winter (December, January and February); MAM, spring (March, April and May); JJA, summer (June, July and August); SON, autumn (September, October and November). The variables represented in each of the first three functions are marked with asterisks. Discriminant function

Precipitation (DJF) Precipitation (MAM) Precipitation (JJA) Precipitation (SON) Annual precipitation Temperature (DJF) Temperature (MAM) Temperature (JJA) Temperature (SON) Mean annual temperature Water balance (DJF) Water balance (MAM) Water balance (JJA) Water balance (SON) Annual water balance Latitude Elevation Evergreen tree species Deciduous tree species

1

2

3

−0.258 0.211 0.557* 0.172 0.281 −0.071 −0.138 −0.304 −0.165 −0.156 −0.274 0.370 0.615* 0.323 0.395 −0.108 −0.187 0.471 −0.471

0.174 0.065 −0.137 0.111 0.028 −0.007 0.016 0.094 0.011 0.021 0.080 0.031 −0.190 0.057 −0.440 −0.134 0.051 0.490* −0.490*

0.373* 0.491* 0.402 0.778* 0.606* 0.372 0.279 0.055 0.283 0.291 0.334 0.392* 0.410 0.769* −0.166 −0.068 0.078 −0.297 0.297

functions 2 and 3. Component 2 (+) and component 3 (−) were characterized by similar centroids in functions 2 and 3. This pattern was also observed for component 2 (−) and component 3 (+), indicating that component 2 (+) and component 3 (−) correspond mostly to evergreen forests. Moreover, component 2 (−) and component 3 (+) were located in more humid sites than those of component 2 (+) and component 3 (−). These latter groups are associated with forests that receive less precipitation at the annual and seasonal scales, particularly during the boreal autumn. They are also located in colder regions, with respect to those of component 2 (−) and component 3 (+). Table 2 summarizes the dominant environmental characteristics controlling the different patterns of response of tree-ring growth to NDVI variations, as detailed in Section 4.2.

5. Discussion and conclusions This study analyzed the relationships between forest growth (tree-ring width indices) and the Normalized Difference Vegetation Index (NDVI) over almost three decades for selected 497 forests across the world. To our knowledge, this is the most comprehensive analysis exploring these relationships at a global scale and covering a wide range of tree species and phenologies.

1.0

1.0

0.8

0.8

0.6

0.6

0.2 0.0

3(-) 2(+)

-0.2

1(-) 4(-)

-0.4

4(+)

0.4

1(+)

FUNCTION 3

FUNCTION 2

2(-)

3(+)

0.4

25

3(+)

2(-)

0.2 0.0

1(+)

1(-)

4(-)

3(-)

-0.2

2(+) -0.4

-0.6

-0.6

-0.8

-0.8

4(+)

-1.0 -1.0 -0.8 -0.6 -0.4 -0.2

0.0

0.2

0.4

0.6

0.8

1.0

FUNCTION 1

-1.0 -1.0 -0.8 -0.6 -0.4 -0.2

0.0

0.2

0.4

0.6

0.8

1.0

FUNCTION 1

1.0 0.8 0.6

4(+)

FUNCTION 3

0.4 0.2

4(-)

0.0

1(+)

-0.2 -0.4

3(+) 2(-)

1(-)

3(-) 2(+)

-0.6 -0.8 -1.0 -1.0 -0.8 -0.6 -0.4 -0.2

0.0

0.2

0.4

0.6

0.8

1.0

FUNCTION 2 Fig. 14. Centroids of the groups obtained through a Principal Components Analysis corresponding to the first three functions of the Predictive Discriminant Analysis (PDA).

26

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Table 2 Dominant environmental characteristics at which the different patterns of response of tree-ring growth to the NDVI are found. Component

Time-scale

Thermal conditions

Hydric conditions

1(+) 1(−) 2(+) 2(−) 3(+) 3(−) 4(+) 4(−)

Short Long Long Short Intermediate Intermediate Short Intermediate

Temperate Temperate Cold Cold Temperate Temperate Temperate Cold

Dry Humid Humid Intermediate Dry Dry Intermediate Humid

5.1. General tree-ring growth vs. NDVI associations We found a general positive relationship between the inter-annual NDVI variability and the annual tree growth in most of the analyzed forests. We also showed that, in almost 67% of the forests, the correlation was positive and statistically significant (p b 0.05). Given the strong relationship between Net Primary Production (NPP) and forest growth, as reported in several studies (e.g., Poulter et al., 2013; Babst et al., 2013; Vicente-Serrano et al., 2015), and also given that a wide range of studies have shown that NDVI is a very good proxy of the forest NPP (e.g., Asrar et al., 1984; Tucker and Sellers, 1986; Coops, 1999), a general positive and significant relationship between NDVI and tree growth can be expected. This finding concurs with previous studies (e.g., Malmstrom et al., 1997; Lopatin et al., 2006; Leavitt et al., 2008; Liang et al., 2009; Lloyd et al., 2010; Macias-Fauria et al., 2012). Nevertheless, the positive and significant relationship locally observed between radial growth and NDVI cannot be considered a rule at the global scale. Other studies have shown no significant relationship between NDVI and tree-ring width (D'Arrigo et al., 2000; Beck et al., 2013). For example, Beck et al. (2013) showed -at four arctic tree-line sites in North America- that the NDVI did not correlate significantly with tree-ring width. Similarly, in mature spruce trees, Beck et al. (2011) demonstrated that changes in the NDVI since 1982 were not consistently correlated with high-frequency changes in tree-ring width. In addition, Berner et al. (2013) documented cases of decreasing tree-ring widths in areas of increasing summer NDVI near a tree-line located in northern Russia; a similar pattern to what was found by Lawrence et al. (2005). Most of these studies based their analysis on the correlations between the tree-ring and the NDVI over specific periods (biweekly or monthly periods). Nevertheless, this approach could not be the optimum to identify a relationship between NDVI and tree growth. This is particularly because some studies have stressed that the annual wood formation is related more to the cumulative NPP recorded over long periods. Huston and Wolverton (2009) stressed that the biomass gain and increase in the long-term stored carbon does not depend on the photosynthetic rates, but rather on how trunks and branches as well as below ground components grow and form wood. This makes the storage of long-lived forms of biomass in the woody parts of trees physically and temporally separated from the leaves carbon gain, which is measured by satellites and quantified in the form of NDVI. Thus, growth and NPP are tightly coupled at annual scales. However, this association may disappear at monthly or daily scales (Granier et al., 2008; Lempereur et al., 2015), a finding confirmed by Zweifel et al. (2010) using automatic dendrometers in a Norway spruce (Picea abies) forest in the Swiss Alps. This behavior could be related to the fact that wood production is the result of accumulating the surplus of synthesized carbohydrates (i.e., it is a secondary carbon sink), and therefore, secondary growth and carbon storage reflect cumulative NPP (Gough et al., 2008). This pattern can be expected because growth and NPP may not be coupled at short temporal scales because wood formation is just one aspect of tree growth and carbon must first be used for primary growth in

order to form shoots, buds, leaves and roots (Stoy et al., 2009). In addition, temporal lags exist during xylogenesis from the expansion to the lignification of wood cells (Cuny et al. 2015). All these processes could explain why we found a general positive and significant correlation in most of the forests analyzed using different time-spans of NDVI accumulation. Here, we must also stress the strong differences found in the tree growth responses to NDVI at the global scale. This feature can be expected given that the annual tree-ring records in each forest respond in a different way to the magnitude, seasonality and accumulation period of the NDVI values. Furthermore, the time period during which plants are actually growing varies dramatically over the Earth, with almost twelve months near the equator, and a month or even less at high latitudes, most elevated regions (Huston and Wolverton, 2009; Cuny et al., 2015) as well as dry sites (Schwinning and Sala, 2004). 5.2. Diversity of the response patterns at the global scale We found eight main patterns of tree-ring response to the NDVI. The first one (component 1 (+)) showed a response to short cumulative NDVI values during the boreal summer. Most studied forests corresponded to this pattern, albeit with no clear geographic structure since this component is characteristic of deciduous forests affected by summer dryness. These forests are impacted by the most negative water balance conditions among all forests analyzed in this study. Different dendrochronological studies focusing on arid and warm forests have shown that xylem production is confined to one season in most temperate and boreal forests (Larson, 1994). Although bimodal xylogenesis (i.e., the production of xylem in two different seasonal periods) has also been described in the Mediterranean and tropical woody species (see, for example, Venugopal and Krishnamurthy, 1987; De Luis et al., 2007; Camarero et al., 2010), xylem production would be constrained by soil water availability as well as photoperiod. These features could make the link between growth and photosynthetic activity restricted more to short periods in which water is available, suggesting that the NDVI is related to growth only during short cumulative periods. Also, a high percentage of forests were represented by a pattern characterized by high correlations with the NDVI values recorded over the previous year (e.g., component 1(−)). These forests correspond mostly to evergreen forests located in humid and temperate areas, such as northwestern portions of North America. This pattern could be explained by the widely known synthesis and storage of carbohydrates by trees, where the tree-activity affects tree-ring width of the following growing season (Kagawa et al., 2005). The lag of growth response to previous NDVI can be very long. Babst et al. (2014a, 2014b) suggested that carbon sequestered after June/July in temperate forests is mostly used for cell-wall thickening processes and/or stored in above- and below-ground nonstructural carbohydrate reserves, which would support next year spring growth (Skomarkova et al., 2006). Accordingly, Richardson et al. (2013) suggested that the prediction of inter-annual variations in wood growth in temperate forests can be improved by accounting for the mobilization of carbon storage pools that were formed several years earlier. This study reported two relevant patterns of response in the forests of component 3. In particular, these forests were characterized by the highest response of the tree-ring growth to the cumulative NDVI over periods of 2–3 months. Nevertheless, the seasonal pattern was quite different between component 3 (+), with a peak of response at the end of the boreal spring, and component 3 (−), whose peak of response occurred at the end of the boreal summer. Both patterns are representative of the dry evergreen forests, although forests of component 3 (+) are found in warmer sites than those of component 3 (−). Temperature differences between these two types of response could explain the observed patterns, since the period of vegetation activity would be longer in forests of warm climate (component 3 (+)) than those of cold

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

climate (component 3 (−)). The latter would be affected by shorter vegetative periods, which would explain why the influence of the NDVI conditions is maximized only during the boreal summer period. Some specific patterns (component 2 (−)) showed a very small response of tree growth to the NDVI, with significant effects being constrained to the end of the vegetation activity period. This pattern was mainly identified in cold biomes from mountain areas (e.g., the Alps) and in high latitudes of North America. In these areas, snow coverage, low air temperature, and short photoperiod combined together constrain the periods of vegetation activity and subsequently limit the duration of tree-ring formation to summer months (Vaganov et al., 1999; Bergeron et al., 2007). In boreal forests, Kaufmann et al. (2004) suggested that NDVI proxies the physiological status of trees, demonstrating that the status of the canopy during summer has a larger effect on tree vigor than the duration of the canopy leaves. This aspect could explain the punctual response of growth to the NDVI that characterizes these forests. 5.3. Summary and concluding remarks We found a very high spatial and temporal diversity in the responses of forest growth to NDVI at the global scale. This reinforces the need of calculating cumulative NDVI values to find a coherent and robust relationship with tree-ring records. The physiological and phenological processes that account for the dependency between wood formation, photosynthetic processes behind the NDVI and time-spans and periods in which these processes are recorded are very complex and not adequately approached. For example, Cuny et al. (2012) showed that growing seasons of wood formation and the resulting growth rates differ between trees and species coexisting even in the same site. Similar patterns were obtained by Camarero et al. (2010) for species located in the same geographical domain, which were described as bimodal and unimodal xylogenesis patterns and explained by local climatic conditions. In addition, the use of the carbon sequestered by trees may vary significantly between seasons, impacting the response of the tree-ring growth to the NDVI. This behavior would really make it difficult to properly unveil the mechanisms that link NDVI variability to radial growth. Moreover, climate constraints (mostly low temperature and humidity) and photoperiod may interfere the relationship between photosynthetic activity and forest growth, since xylogenesis can be recorded in very short periods of high growth rates, in response to light availability or temperature conditions (Mäkinen et al., 2003; Kirdyanov et al., 2003; Rossi et al., 2006). Thus, the different climate conditions demanded for wood formation and leaf activity may also uncouple NDVI-growth relationship, as higher temperatures are required for growth than for photosynthesis (Körner, 1998; Suni et al., 2003; Rossi et al., 2011). Although the present study is constrained by the temporal availability of data (around 30-years), the spatial distribution of the available samples - given the scarce number of samples in the tropical forestsand the availability of the remote-sensing information that is currently coarse over space and often noisy, our findings can provide a general methodological and conceptual frame for detailed regional studies. In the same context, while some uncertainties might be introduced in this study, given the chosen tree-ring database and the coarse resolution of the NDVI data, we still believe that this study provides the first global assessment of the response of tree-ring growth to the NDVI, introducing useful and novel results with potential implications for understanding forest growth under different environmental conditions. From a practical perspective, to better link the findings of this work with different ecological applications at detailed spatial scales (e.g., forest management, developing transfer functions to predict the tree-ring growth, etc.), more detailed analysis would be needed using satellite imagery of higher spatial and temporal resolution, in addition to a dense and representative number of collected samples at the sites of interest. Although in this study we have showed that the robustness of the results obtained from low-resolution GIMMS NDVI is high since they are

27

comparable with those obtained in most of the forest with high resolution MODIS data (1 km), it is expected that higher resolution spacebased products (e.g., MODIS) can allow a better investigation of the relationship between tree-ring growth and NDVI in different environments and forest types, particularly when longer time series become available.

Acknowledgements This work was supported by the research projects PCIN-2015-220, CGL2014-52135-C03-01, CGL2015-69186-C2-1-R Red de variabilidad y cambio climático RECLIM (CGL2014-517221-REDT), Red de excelencia ecometas (CGL2014-53840-REDT) financed by the Spanish Commission of Science and Technology and FEDER funds, “IMDROFLOOD-Improved Drought and Flood Early Warning, Forecasting and Mitigation using real-time hydroclimatic indicators” supported under the Water Joint programme Initiative “Water Challenges for a Changing World”, Water Works 2014 Cofunded call and “LIFE12 ENV/ES/000536-Demonstration and validation of innovative methodology for regional climate change adaptation in the Mediterranean area (LIFE MEDACC)” financed by the LIFE programme of the European Commission. A. Gazol is supported by a Postdoctoral grant from MINECO (Contrato Formación Postdoctoral MINECO - FPDI 2013-16600, FEDER funds).

References Allen, C.D., et al., 2010. A global overview of drought and heat-induced tree mortality reveals emerging climate change risks for forests. For. Ecol. Manag. 259, 660–684. Asrar, G., Fuchs, M., Kanemasu, E.T., Hatfield, J.L., 1984. Estimating absorbed photosynthetic radiation and leaf area index from spectral reflectance in wheat. Agron. J. 76, 300–312. Babst, F., Poulter, B., Trouet, V., et al., 2013. Site- and species-specific responses of forest growth to climate across the European continent. Glob. Ecol. Biogeogr. 22, 706–717. Babst, F., Bouriaud, O., Alexander, R., Trouet, V., Frank, D., 2014a. Toward consistent measurements of carbon accumulation: a multi-site assessment of biomass and basal area increment across Europe. Dendrochronologia 32, 153–161. Babst, F., Bouriaud, O., Papale, D., et al., 2014b. Above-ground woody carbon sequestration measured from tree rings is coherent with net ecosystem productivity at five eddycovariance sites. New Phytol. 201, 1289–1303. Barber, V.A., Juday, G.P., Finney, B.P., 2008. Reduced growth of Alaskan white spruce in the twentieth century from temperature-induced drought stress. Glob. Planet. Chang. 60, 289–305. Baret, F., Guyot, G., 1991. Potentials and limits of vegetation indices for LAI and APAR assessment. Remote Sens. Environ. 35, 161–173. Barry, R.G., Carleton, A.M., 2001. Synoptic and Dynamic Climatology. Routledge, London. Beck, P.S.A., et al., 2011. Satellite observations of high northern latitude vegetation productivity changes between 1982 and 2008: ecological variability and regional differences. Environ. Res. Lett. 6, 045501. Beck, P.S.A., Andreu-Hayles, L., D'Arrigo, E., et al., 2013. A large-scale coherent signal of canopy status in maximum latewood density of tree rings at arctic treeline in North America. Glob. Planet. Chang. 100, 109–118. Bergeron, O., et al., 2007. Comparison of carbon dioxide fluxes over three boreal black spruce forests in Canada. Glob. Chang. Biol. 13, 89–107. Berner, L.T., Beck, P.S.A., Bunn, A.G., Goetz, S.J., 2013. Plant response to climate change along the forest-tundra ecotone in northeastern Siberia. Glob. Chang. Biol. 19, 3449–3462. Bigler, C., Bräker, O.U., Bugmann, H., Dobbertin, M., Rigling, A., 2006. Drought as an inciting mortality factor in scots pine stands of the Valais, Switzerland. Ecosystems 9, 330–343. Bontemps, S., Defourny, P., Radoux, J., Van Bogaert, E., Lamarche, C., Achard, F., Mayaux, P., Boettcher, M., Brockmann, C., Kirches, G., Zülkhe, M., Kalogirou, V., Arino, O., 2013. Consistent Global Land Cover Maps for Climate Modeling Communities: Current Achievements of the ESA's Land Cover CCIESA Living Planet Symposium 9–13 September 2013. Edinburgh, United Kingdom. Büntgen, U., Trouet, V., Frank, D., et al., 2010. Tree-ring indicators of German summer drought over the last millennium. Quat. Sci. Rev. 29, 1005–1016. Camarero, J.J., Albuixech, J., López-Lozano, R., Casterad, M.A., 2010. An increase in canopy cover leads to masting in Quercus ilex. Trees 24, 909–918. Camarero, J.J., Franquesa, M., Sangüesa-Barreda, G., 2015. Timing of drought triggers distinct growth responses in holm oak: implications to predict warming-induced forest defoliation and growth decline. Forests 6, 1576–1597. Carlson, T.N., Ripley, D.A., 1997. On the relation between NDVI, fractional vegetation cover, and leaf area index. Remote Sens. Environ. 62, 241–252. Cihlar, J., StLaurent, L., Dyer, J.A., 1991. The relation between normalized difference vegetation index and ecological variables. Remote Sens. Environ. 35, 279–298. Cook, E.R., Kairiukstis, L.A., 1990. Methods of dendrochronology: applications in the environmental sciences.

28

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29

Cook, E.R., Krusic, P.J., 2005. ARSTAN, A Tree-ring Standardization Program Based on Detrending and Autoregressive Time Series Modeling, with Interactive Graphics. Lamont-Doherty Earth Observatory of Columbia University, Palisades, New York, U.S.A, Tree-Ring Laboratoryhttp://www.ldeo.columbia.edu/tree-ring-laboratory/ resources/software. Coops, N., 1999. Improvement in predicting stand growth of Pinus radiafa (D. Don) across landscapes using NOAA AVHRR and Landsat MSS imagery combined with a forest growth process model (3-PGS). Photogramm. Eng. Remote Sens. 65, 1149–1156. Crowther, T.W., Glick, H.B., Covey, K.R., et al., 2015. Mapping tree density at a global scale. Nature 525, 201–205. Cuny, H.E., Rathgeber, C.B.K., Lebourgeois, F., Fortin, M., Fournier, M., 2012. Life strategies in intra-annual dynamics of wood formation: example of three conifer species in a temperate forest in north-east France. Tree Physiol. 32, 612–625. Cuny, H.E., Rathgeber, C.B.K., Frank, D., et al., 2015. Woody biomass production lags stemgirth increase by over one month in coniferous forests. Nature Plants 1, 1–6. D'Arrigo, R.D., Malmstrom, C.M., Jacoby, G.C., Los, S.O., Bunker, D.E., 2000. Correlation between maximum latewood density of annual tree rings and NDVI based estimates of forest productivity. Int. J. Remote Sens. 21, 2329–2336. De Luis, M., Gričar, J., Čufar, K., Raventós, J., 2007. Seasonal dynamics of wood formation in Pinus halepensis from dryand semi-arid ecosystems in Spain. IAWA J. 28, 389–404. Dobbertin, M., 2005. Tree growth as indicator of tree vitality and of tree reaction to environmental stress: a review. Eur. J. For. Res. 124, 319–333. Duncan, J., Stow, D., Fanklin, J., Hope, J., 1993. Assessing the relationship between spectral vegetation indices and shrub cover in the Jornada Basin, New Mexico. Int. J. Remote Sens. 14, 3395–3416. Fang, J., Chen, A., Peng, C., Zhao, S., Ci, L., 2001. Changes in forest biomass carbon storage in China between 1949 and 1998. Science 292, 2320–2322. Fridman, J., Walheim, M., 2000. Amount, structure, and dynamics of dead wood on managed forestland in Sweden. For. Ecol. Manag. 131, 23–36. Fritts, H.C., 2001. Tree Rings and Climate. Blackburn Press 567 pp. Gillies, R.R., Carlson, T.N., Cui, J., Kustas, W.P., Humes, K.S., 1997. A verification of the triangle method for obtaining surface soil water content and energy fluxes from remote measurements of the normalized difference vegetation index (NDVI) and surface radiant temperature. Int. J. Remote Sens. 18, 3145–3166. Gough, C.M., Vogel, C.S., Schmid, H.P., Su, H.B., Curtis, P.S., 2008. Multi-year convergence of biometric and meteorological estimates of forest carbon storage. Agric. For. Meteorol. 148, 158–170. Gouveia, C.M., Páscoa, P., Russo, A., Trigo, R.M., 2016. Land degradation trend assessment over iberia during 1982-2012. Cuadernos de Investigacion Geografica 42 (1), 89–112. Goward, S.N., Dye, D.G., 1987. Evaluating North American net primary productivity with satellite observations. Adv. Space Res. 7, 165–174. Granier, A., Bréda, N., Longdoz, B., Gross, P., Ngao, J., 2008. Ten years of fluxes and stand growth in a young beech forest at Hesse, North-eastern France. Ann. For. Sci. 65, 704. Grissino-Mayer, H.D., Fritts, H.C., 1997. The international tree-ring data bank: an enhanced global database serving the global scientific community. The Holocene 7, 235–238. Gutman, G.G., 1991. Vegetation indices from AVHRR. An update and future prospects. Remote Sens. Environ. 35, 121–136. Hair, J.F., Anderson, R.E., Tatham, R.L., Black, W.C., 1998. Multivariate data analysis. Prentice Hall, New York, NY. Hall, F.G., Botkin, D.B., Strebel, D.E., Woods, K.D., Goetz, S.J., 1991. Large-scale patterns of forest succession as determined by remote sensing. Ecology 72, 628–640. Harris, I., Jones, P.D., Osborn, T.J., Lister, D.H., 2014. Updated high-resolution grids of monthly climatic observations - the CRU TS3.10 dataset. Int. J. Climatol. 34, 623–642. Hartmann, H., Ziegler, W., Kolle, O., Trumbore, S., 2013. Thirst beats hunger - declining hydration during drought prevents carbon starvation in Norway spruce saplings. New Phytol. 200, 340–349. Hasenauer, H., Petritsch, R., Zhao, M., Boisvenue, C., Running, S.W., 2012. Reconciling satellite with ground data to estimate forest productivity at national scales. For. Ecol. Manag. 276, 196–208. Hirota, M., Holmgren, M., Van Nes, E.H., Scheffer, M., 2011. Global resilience of tropical forest and savanna to critical transitions. Science 334, 232–235. Huberty, C.J., 1994. Applied Discriminant Analysis. Wiley, New York, NY. Huston, M.A., Wolverton, S., 2009. The global distribution of net primary production: resolving the paradox. Ecol. Monogr. 79, 343–377. Jenkins, J.C., Chojnacky, D.C., Heath, L.S., Birdsey, R.A., 2003. National-scale biomass estimators for United States tree species. For. Sci. 49, 12–35. Kagawa, A., Sugimoto, A., Yamashita, K., Abe, H., 2005. Temporalphotosynthetic carbon isotope signatures revealed in a tree ring through 13CO2 pulse-labelling. Plant Cell Environ. 28, 906–915. Kaufmann, R.K., D'Arrigo, R.D., Laskowski, C., Myneni, R.B., Zhou, L.M., Davi, N.K., 2004. The effect of growing season and summer greenness on northern forests. Geophys. Res. Lett. 31, L09205. http://dx.doi.org/10.1029/2004GL019608. Kaufmann, R.K., D'Arrigo, R.D., Paletta, L.F., Tian, H.Q., Jolly, W.M., Myneni, R.B., 2008. Identifying climatic controls on ring width: the timing of correlations between tree rings and NDVI. Earth Interactions 12, 14. Ketterings, Q.M., Coe, R., Van Noordwijk, M., Ambagau, Y., Palm, C.A., 2001. Reducing uncertainty in the use of allometric biomass equations for predicting above-ground tree biomass in mixed secondary forests. For. Ecol. Manag. 146, 199–209. Kirdyanov, A., Hughes, M., Vaganov, E., Schweingruber, F., Silkin, P., 2003. The importance of early summer temperature and date of snow melt for tree growth in the Siberian subarctic. Trees - Structure and Function 17, 61–69. Körner, C., 1998. A re-assessment of high elevation treeline positions and their explanation. Oecologia 115, 445–459. Larson, P.R., 1994. The Vascular Cambium: Development and Structure. Springer-Verlag, Heidelberg, Germany.

Lawrence, G.B., Lapenis, A.G., Berggren, D., et al., 2005. Climate dependency of tree growth suppressed by acid deposition effects on soils in Northwest Russia. Environ. Sci. Technol. 39, 2004–2010. Leavitt, S.W., Chase, T.N., Rajagopalan, B., Lee, E., Lawrence, P.J., 2008. Southwestern U.S. tree-ring carbon isotope indices as a possible proxy for reconstruction of greenness of vegetation. Geophys. Res. Lett. 35 (L12704). http://dx.doi.org/10.1029/ 2008GL033894. Lempereur, M., Martin-StPaul, N.K., Damesin, C., Joffre, R., Ourcival, J.M., Rocheteau, A., Rambal, S., 2015. Growth duration is a better predictor of stem increment than carbon supply in a Mediterranean oak forest: implications for assessing forest productivity under climate change. New Phytol. 207, 579–590. Liang, E., Eckstein, D., Liu, H., et al., 2009. Assessing the recent grassland greening trend in a long-term context basedon tree-ring analysis: a case study in North China. Ecol. Indic. 9, 1280–1283. Liu, H., Park Williams, A., Allen, C.D., et al., 2013. Rapid warming accelerates tree growth decline in semi-arid forests of Inner Asia. Glob. Chang. Biol. 19, 2500–2510. Lloyd, A.H., Fastie, C.L., 2002. Spatial and temporal variability in the growth and climate response of treeline trees in Alaska. Clim. Chang. 52, 481–509. Lloyd, A.H., Bunn, A.G., Berner, L., 2010. A latitudinal gradient in tree growth response toclimate warming in the Siberian taiga. Glob. Chang. Biol. 17, 1935–1945. Lopatin, E., Kolstrom, T., Spiecker, H., 2006. Determination of forest growth trends in Komi Republic (northwestern Russia): combination of tree-ring analysis andremote sensing data. Boreal Environ. Res. 11, 341–353. Lu, D., 2006. The potential and challenge of remote sensing-based biomass estimation. Int. J. Remote Sens. 27, 1297–1328. Macias-Fauria, M., Forbes, B.C., Zetterberg, P., Kumpula, T., 2012. Eurasian Arctic greening reveals teleconnections and the potential for structurally novel ecosystems. Nat. Clim. Chang. 2, 613–618. Mäkinen, H., Nöjd, P., Saranpää, P., 2003. Seasonal changes in stem radius and production of new tracheids in Norway spruce. Tree Physiol. 23, 959–968. Malmstrom, C.M., Thompson, M.V., Juday, G., Los, S.O., Randerson, J.T., Field, C.B., 1997. Interannual variation in global-scale net primary production: testing modelestimates. Glob. Biogeochem. Cycles 11, 367–392. Martínez-Vilalta, J., Piñol, J., 2002. Drought-induced mortality and hydraulic architecture in pine populations of the NE Iberian Peninsula. For. Ecol. Manag. 161, 247–256. Myneni, R.B., Hall, F.G., Sellers, P.J., Marshak, A.L., 1995. The meaning of spectralvegetation indices. IEEE Trans. Geosci. Remote Sens. 33, 481–486. Pasho, E., Camarero, J.J., de Luis, M., Vicente-Serrano, S.M., 2011. Impacts of drought at different time scales on forest growth across a wide climatic gradient in north-eastern Spain. Agric. For. Meteorol. 151, 1800–1811. Peterson, D.W., Peterson, D.L., 2001. Mountain hemlock growth responds to climatic variability at annual and decadal time scales. Ecology 82, 3330–3345. Phdersen, B.S., 1998. The role of stress in the mortality of Midwestern oaks as indicated by growth prior to death. Ecology 9, 79–93. Pinzon, J.E., Tucker, C.J., 2014. A non-stationary 1981–2012 AVHRR NDVI3g time series. Remote Sens. 6, 6929–6960. Poulter, B., Pederson, N., Liu, H., et al., 2013. Recent trends in Inner Asian forest dynamics to temperature and precipitation indicate high sensitivity to climate change. Agric. For. Meteorol. 178–179, 31–45. Richardson, A.D., Carbone, M.S., Keenan, T.F., et al., 2013. Seasonal dynamics and age of stemwood nonstructural carbohydrates in temperate forest trees. New Phytol. 197, 850–861. Richman, M.B., 1986. Rotation of principal components. J. Climatol. 6, 293–335. Rossi, S., Deslauriers, A., Anfodillo, T., et al., 2006. Conifers in cold environments synchronize maximum growth rate of tree-ring formation with day length. New Phytol. 170, 301–310. Rossi, S., Morin, H., Deslauriers, A., Plourde, P.Y., 2011. Predicting xylem phenology in black spruce under climate warming. Glob. Chang. Biol. 17, 614–625. Rouse, J.W., Haas, R.H., Shell, J.A., Deering, D.W., 1974. Monitoring vegetation systems in the Great Plains with ERTS-1. In: SC, F., EP, M., MA, B. (Eds.), Third ERTS-1 Symposium. NASA, Washington, pp. 309–317. Running, S.W., Nemani, R.N., Heinsch, F.A., et al., 2004. A continuous satellite-derived measure of global terrestrial primary production. Bioscience 54, 547–560. Sarris, D., Christodoulakis, D., Körner, C., 2007. Recent decline in precipitation and tree growth in the eastern Mediterranean. Glob. Chang. Biol. 13, 1187–1200. Schwinning, S., Sala, O.E., 2004. Hierarchy of responses to resource pulses in arid and semi-arid ecosystems. Oecologia 141, 211–220. Shabanov, N.V., Huang, D., Yang, W., et al., 2005. Analysis and optimization of the MODIS leaf area index algorithm retrievals over broadleaf forests. IEEE Trans. Geosci. Remote Sens. 43, 1855–1865. Shestakova, T.A., Gutiérrez, E., Kirdyanov, A.V., et al., 2016. Forests synchronize their growth in contrasting Eurasian regions in response to climate warming. Proc. Natl. Acad. Sci. U. S. A. 113, 662–667. Skomarkova, M.V., Vaganov, E.A., Mund, M., et al., 2006. Inter-annual and seasonal variability of radial growth, wood density and carbon isotope ratios in tree rings of beech (Fagus sylvatica) growing in Germany and Italy. Trees - Structure and Function 20, 571–586. Song, C., Woodcook, C.E., 2003. Monitoring forest succession with multitemporal landsat images: factors of uncertainty. IEEE Trans. Geosci. Remote Sens. 41, 2557–2567. Stokes, M.A., Smiley, T.L., 1968. An Introduction to Tree-Ring Dating. University of Chicago Press 1968 - Science - 73 pp. Stoy, P.C., Richardson, A.D., Baldocchi, D.D., 2009. Biosphere-atmosphere exchange of CO2 in relation to climate: a cross-biome analysis across multiple time scales. Biogeosciences 6, 2297–2312. Suni, T., Rinne, J., Reissell, A., et al., 2003. Long-term measurements of surface fluxes above a Scots pine forest in Hyytiälä, southern Finland, 1996–2001. Boreal Environ. Res. 8, 287–301.

S.M. Vicente-Serrano et al. / Remote Sensing of Environment 187 (2016) 14–29 Takao, G., Priyadi, H., Ikbal Nursal, W., 2010. The Operational Role of Remote Sensing in Forest and Landscape Management: Focus Group Discussion Proceedings. Center for International Forestry Research (CIFOR), Bogor, Indonesia 97 pp. Tardif, J., Camarero, J.J., Ribas, M., Gutiérrez, E., 2003. Spatiotemporal variability in tree growth in the Central Pyrenees: climatic and site influences. Ecol. Monogr. 73, 241–257. Timothy, D., Onisimo, M., Cletah, S., Adelabu, S., Tsitsi, B., 2016. Remote sensing of above ground forest biomass: a review. Trop. Ecol. 57, 125–132. Tucker, C.J., 1979. Red and photographic infrared linear combinations for monitoringvegetation. Remote Sens. Environ. 8, 127–150. Tucker, C.J., Sellers, P.J., 1986. Satellite remote sensing of primary production. Int. J.Remote Sens. 7, 1395–1416. Tucker, C.J., Pinzon, J.E., Brown, M.E., Slayback, D.A., Pak, E.W., Mahoney, R., Vermote, E.F., El Saleous, N., 2005. An extended AVHRR 8-km NDVI dataset compatible with MODIS andSPOT vegetation NDVI data. Int. J. Remote Sens. 26, 4485–4498. Vaganov, E.A., Hughes, M.K., Kirdyanov, A.V., Schweingruber, F.H., Silkin, P.P., 1999. Influence of snowfall and melt timing on tree growth in subarctic Eurasia. Nature 400 (6740), 149–151. Van Breugel, M., Ransijn, J., Craven, D., Bongers, F., Hall, J.S., 2011. Estimating carbon stock in secondary forests: decisions and uncertainties associated with allometric biomass models. For. Ecol. Manag. 262, 1648–1657.

29

Venugopal, N., Krishnamurthy, K.V., 1987. Seasonal production of secondary phloem in the twigs of certain tropical timber trees. Ann. Bot. 60, 61–67. Vicente-Serrano, S.M., 2007. Evaluating the impact of drought using remote sensing in a Mediterranean, semi-arid region. Nat. Hazards 40, 173–208. Vicente-Serrano, S.M., Camarero, J.J., Zabalza, J., et al., 2015. Evapotranspiration deficit controls net primary production and growth of silver fir: implications for CircumMediterranean forests under forecasted warmer and drier conditions. Agric. For. Meteorol. 206, 45–54. Williams, A.P., Allen, C.D., Millar, C.I., et al., 2010. Forest responses to increasing aridity and warmth in the southwestern United States. Proc. Natl. Acad. Sci. U. S. A. 107, 21289–21294. Williams, A.P., Allen, C.D., Macalady, A.K., et al., 2013. Temperature as a potent driver of regional forest drought stress and tree mortality. Nat. Clim. Chang. 3, 292–297. Wylie, B.K., Meyer, D.J., Tieszen, L.L., Mannel, S., 2002. Satellite mapping of surface biophysical parameters at the biome scale over the North American grasslands a case study. Remote Sens. Environ. 79, 266–278. Zweifel, R., Eugster, W., Etzold, S., et al., 2010. Link between continuous stem radius changes and net ecosystem productivity of a subalpineNorway spruce forest in the Swiss Alps. The New Phytologist 187, 819–830.