Silver Fir Defoliation Likelihood Is Related to Negative Growth Trends and High Warming Sensitivity at Their Southernmost Distribution Limit

Changes in radial growth have been used to estimate tree decline probability since they may indicate tree responses to long- and short-term stressors. We used visual assessments of crown defoliation, an indicator of decline, and retrospective tree-ring analyses to determine whether climate-growth sensitivity and tree growth rates may be used as predictors of tree die-o ﬀ probability in Abies alba (silver ﬁr) at the Spanish Pyrenees. We used climatic data to calculate standardized temperature and precipitation data and drought indexes. Basal area increment was measured for declining (defoliation > 50%) and nondeclining (defoliation < 50%) silver ﬁrs in stands with contrasting defoliation. Logistic regressions were applied to predict tree die-o ﬀ . Since the early 1980s, a synchronised reduction in basal area increment was observed in declining trees. The basal area increment trend correctly classiﬁed 64% of declining trees and 94% of nondeclining trees. The growth sensitivity to water deﬁcit, temperature, and a drought index also signiﬁcantly predicted silver ﬁr decline, but providing underestimated predictions. Our ﬁndings underscore the idea that long-term climatic warming seems to be a major driver of growth decline in silver ﬁr. Ongoing growth reduction and enhanced mortality may promote vegetation shifts in declining Pyrenean A. alba forests.


Introduction
Declining trends of tree radial growth are considered as reliable indicators of long-term stress and may be an additional risk factor for drought-induced mortality [1][2][3]. Individual tree decline and death often occur as a result of the combined effects of different long-and short-term stressors [4,5]. Radial growth, used here as a proxy of whole plant carbon gain [6], may be used to identify those trees with the highest probabilities of death [3,7]. Indeed, several studies have shown that mortality rates are to some extent inversely related to radial growth trends [8][9][10][11]. However, the potential use of both growth trends and climate sensitivity, as surrogates for evaluating the vulnerability of tree species to climate change, has received less attention [3,[12][13][14].
Dendrochronological assessments of changing trends of radial growth may be useful to understand decline processes [15,16]. Usually, reduced wood formation occurs prior to visual symptoms of decline such as crown defoliation. Thus, dendrochronology may be useful to forecast the impending decline of particular trees and forests [3,16]. Moreover, tree decline may result of the combined effects of several stressors acting at different time scales on growth. Therefore, including information from both long and short-term stressors should improve our knowledge of decline mechanisms taking also into account specific traits of tree species (e.g., drought tolerance).
To evaluate the contributing effects on tree decline of tree growth trend (i.e., contrasting growth patterns obtained in trees subjected to the same regional climate) and climate responsiveness, we used crown defoliation, an indicator of tree decline, and retrospective tree-ring analyses. We focused on declining silver fir (Abies alba Mill.) stands from  the Spanish Pyrenees [17,18]. We hypothesized that silver fir decline is a combined response to climatic stressors, acting at short-and long-term scales, and changes in growth trends, respectively. Thus, an assessment of the effects of both types of stressors on growth should provide improved models of tree decline. We also hypothesized that declining growth trends reduce the capacity of trees to adapt to drier climatic conditions. Our specific aims were to (i) quantify growth patterns in declining and nondeclining A. alba trees, (ii) provide predictors of silver fir decline at the individual level based on growth patterns and growth sensitivity to climate variability.

Study Species and Field
Sampling. The Pyrenees constitute a transitional area between more humid conditions in their northern margin and drier conditions southwards where Mediterranean vegetation becomes dominant [19]. The studied A. alba populations are located in the Aragón Pyrenees, northeastern Spain (Figure 1), where silver fir stands are usually found at humid sites on north-facing slopes forming pure or mixed stands with Fagus sylvatica L. or Pinus sylvestris L. Most studied stands are located on marls and limestones, which generate basic soils, or on moraine deposits with rocky but deep soils. The most used method of timber harvesting in the study area was diameter limit cutting, which mostly affected fast-growing and big trees. According to historical data, logging intensity during the 20th century in the Pyrenees was greatest in the 1950s but no data are available on how widespread was in this region [20].
Since the 1980s silver fir growth decline in the Aragón Pyrenees was characterized by high defoliation levels and increased mortality rates (Figure 2), mainly at low-elevation sites [21]. Since we were interested in comparing A. alba stands with different decline symptoms, we performed an extensive field survey, visiting at least one site in all 10 km 2 grids where silver fir formed forests across the study area. Based on this extensive field survey, we selected representative stands within the coordinates 42 • 42 -42 • 45 N, 0 • 38 -0 • 52 W; 1000-1300 m a.s.l., with about 30% of declining trees (see below the definition of declining and nondeclining trees; detailed study site characteristics can be found in [18]). A total of 62 dominant trees (28 declining and 34 nondeclining silver firs; ∼40 cm diameter from the base, ∼20 m height) were selected for sampling within a 500 m long and 20 m wide transect randomly located within the stands.

Climate Data.
We used local climatic records to study the spatiotemporal variation of climatic conditions in the study area (detailed characteristics of climate data can be found in [23]). To estimate the missing data for each meteorological station and to obtain a regional climatic record we used the program MET from the Dendrochronology Program Library [24]. The annual water budget (named as drought index thereafter; P-ETP) was obtained from the sum of the differences between monthly data of precipitation (P) and potential evapotranspiration (ETP), and the last variable was estimated by a modified version of the Thornthwaite method [25].
To obtain the regional climatic mean local climate data were standardized to give each station the same relative weight. Mean temperature and P-ETP were standardized by subtracting the mean and dividing by the standard deviation: where Z i is the x i score distance from the x average (x) in standard deviation units (σ).
To determine the severity and the statistical rarity of extreme drought events, precipitation data was standardized by calculating the standardized precipitation index (SPI) for each meteorological station. The spatial and temporal dimensions of drought events create problems in generating a drought index because not only must an anomaly be normalized with respect to location, but the anomaly must also be normalized in time if it is to produce a meaningful estimate of drought. The SPI accomplishes both. The SPI is normalized to a station location because it accounts for the frequency distribution of precipitation as well as the accompanying variation at the station. Additionally, the SPI is normalized in time because it can be computed at any number of time scales, depending upon the impacts of interest to the analyst. Conceptually, the SPI represents a z-score, or the number of standard deviations above or below that an event is from the mean, as it is expressed in (1). However, the SPI performs a preadjustment to this standard formulation due to precipitation is typically positively skewed. To adjust for this mathematical rainfall feature, the precipitation data is transformed to a more normal distribution by applying the gamma function [26][27][28]. The gamma distribution is defined by its frequency or probability density function: where α > 0, is a shape parameter; β > 0, is a scale parameter; x > 0, is the precipitation amount; Γ(α) is the gamma function: Computation of the SPI involves fitting a gamma probability density function to a given frequency distribution of precipitation totals for a station. The alpha and beta parameters of the gamma probability density function are estimated for each station, for each time scale of interest (seasonal, annual, etc.), and for each month of the year. The maximum likelihood solutions are used to optimally estimate α and β as follows: where where n = number of precipitation observations. The resulting parameters are then used to find the cumulative probability of an observed precipitation event for the given month and time scale for the station in question. The cumulative probability is given by: Letting t = x/ β, this equation becomes the incomplete gamma function: Since the gamma function is undefined for x = 0 and a precipitation distribution may contain zeros, the cumulative probability becomes: where q is the probability of a zero. If m is the number of zeros in a precipitation time series, q can be estimated by m/n. The cumulative probability, H(x), is then transformed to the standard normal random variable Z with mean zero and variance of one (see (1)) as follows: where

Dendrochronological Methods.
Dendrochronological sampling was performed in the 62 trees selected following standard methods [29]. Two or three cores were taken from each tree at breast height (1.3 m) using an increment borer. The wood samples were air-dried and polished with a series of successively finer sand-paper grits until tree rings were clearly visible. Then, the wood samples were visually crossdated. Tree rings were measured to the nearest 0.001 mm using a binocular scope and a LINTAB measuring device (Rinntech, Heidelberg, Germany). The visual cross-dating of the tree-rings was checked using the program COFECHA [30]. The trend due to the geometrical constraint of adding a volume of wood to a stem of increasing radius was corrected by converting tree-ring widths into basal area increments (BAI) [31]. We calculated mean BAI values and trends separately for declining (defoliation > 50%) and nondeclining (defoliation < 50%) trees. Trends values were calculated as the slope of the time series, estimated by least squares regression. Growth reduction percentage was calculated as GRP = 100 × ((mean BAI of the last decade included in the compute) − (mean BAI of the first decade Table 1: Basal area increment (BAI) mean values and trends measured in declining (defoliation > 50%) and non-declining (defoliation < 50%) trees in the studied stands. Trends values are calculated as the slope of the time series, estimated by least squares regression. Growth reduction percentage is calculated as GRP = 100 × ((mean BAI of the last decade included in the compute) − (mean BAI of the first decade included in the compute))/(mean BAI of the first decade included in the compute). The coefficient of variation (CV) is defined as the ratio of the standard deviation to the mean multiplied by 100. Different letters indicate significant differences based on for one-way ANOVAs. included in the compute))/(mean BAI of the first decade included in the compute). The coefficient of variation (CV) was defined as the ratio of the standard deviation to the mean multiplied by 100 (Table 1). BAI data were also standardized, following the same procedure described above for climate data [31]. Growth sensitivity to climate was estimated, for each individual tree, as the determination coefficient between standardized basal area increment and annual climatic variables (mean annual temperature, annual SPI, and annual water budget).

Logistic Regression Models of Silver Fir Decline.
We applied logistic regression to predict the decline probability P(Y i,t = 1 | X i,t ) of tree i at a given time t [32]: where Y i,t = 1 indicates that tree i is declining at time t (correspondingly, Y i,t = 0 indicates that tree i is nondeclining at time t). The matrix X i,t contains the independent variables of tree i at time t and β is a vector containing the regression coefficients. The predictors used corresponded to three groups of variables related to growth level, growth trends, and growth sensitivity to climate. The assessed variables were the mean tree BAI for the periods 1960-1999, 1960-1979, and 1970-1999; the BAI trend of each tree computed as the time/BAI slope from least squares linear regression, for the same periods mentioned before; the coefficient of variation of BAI for each tree also calculated for the same periods as above; the growth sensitivity to changes in temperature, precipitation, and drought (P-ETP). Logistic regressions based on generalized linear models for binomial data were fitted by maximum log-likelihood estimation using the R software [33]. The P value and the odds ratio were used during the model selection procedure to assess the goodnessof-fit of the models [34,35].

Results
Standardized mean annual temperature yielded systematically above average values after 1981, with extreme warm years in the following periods 1937-1938, 1986-1990, and 1997-1998    3: Regional standardized mean temperature (T), standardized precipitation index (SPI), standardized total precipitation minus evapotranspiration (P-ETP; difference between precipitation-P-and potential evapotranspiration-ETP), and standardized basal area increment chronologies of declining and nondeclining silver firs. Pearson correlation coefficients (r) for the relationships between standardized basal area increment and the annual climatic variables are also displayed with their significance levels ( * P < 0.05; * * P < 0.01).
The standardized drought index (P-ETP) yielded extreme drought events with pronounced water deficit in 1933, 1937, 1949, and 1989, while years with high water availability were the same that those estimated using the SPI. The standardized BAI were significantly correlated to standardized climate series for both nondeclining and declining trees (Figure 3). The highest correlation coefficient was Table 2: Logistic regression models for declining (n = 28) and non-declining (n = 34) silver firs. The effects of growth characteristics (mean, trends, and coefficient of variation-CV) and growth sensitivity to climate (estimated as the determination coefficient between basal area increment-BAI-and annual climatic variables). The annual water budget (P-ETP) was the differences between precipitation (P) and potential evapotranspiration (ETP). The relative maximum likelihood is calculated as the percentage of maximum likelihood value obtained. Nonsignificant models are indicated by ns.

Model
Maximum obtained between standardized mean annual temperature and BAI (r = −0.41 and r = −0.34 for declining and nondeclining trees, resp.), followed by those obtained between the standardized mean drought index and BAI (r = 0.27 and r = 0.23 for declining and nondeclining trees, resp.). The standardized BAI was also related to the extreme drought events recorded by the SPI series, in the case of declining trees. In the 1950s and late 1970s, higher water availability was related to positive standardized BAI, especially in the case of nondeclining trees. Contrastingly, in 1986 all trees showed a sudden growth decline ( Figure 3). However, only declining trees showed negative standardized BAI between 1995 and 1999. The combination of two or more independent variables increased the probability and the likelihood values of logistic regression models of silver fir decline, but it did not improve the percentage of correctly assigned trees to the declining and non-declining classes (data not shown). The best decline predictor was the BAI trend for the period 1980-1999, followed by growth sensitivity to drought ( Table 2; Figures  4 and 5). The BAI trend calculated for the period 1960-1999, as well as the BAI coefficient of variation estimated for the period 1980-1999, also obtained substantial support. The BAI trend calculated for the period 1980-1999 correctly classified 64% and 94% of declining and non-declining firs, respectively. The sensitivities of growth to changes in drought, temperature, and precipitation resulted in an underestimation of correctly classified declining trees (46%, 27%, and 27%, respectively; Table 2 and Figure 4).

Discussion
Logistic regression modelling showed a relatively high accuracy to predict silver fir decline using the trend in basal area increment as predictor. Logistic regression models based on the slope of the basal area increment of the five years before tree death also showed high predictive ability of competitionmediated mortality of Norway spruce (Picea abies) in the Swiss Alps [32]. In another study [16] decline and crown defoliation in silver firs from Slovenia was modelled, showing that declining trees had a higher probability to die than nondeclining trees.
Defoliated trees followed diminishing growth trends up two decades before presenting external decline symptoms as abundant needle loss, which agrees with other studies [36][37][38], and suggest that responses of trees to climatic stress can be predicted based on characteristics growth patterns [5,32,39]. We were able to detect a significant difference between long-term growth patterns of healthy and declining silver firs. Moreover, in many cases, declining trees showed negative growth trends, and high growth sensitivity to climate. Abrupt declines in growth or strongly negative growth trends may indicate a rapid physiological adaptation to changing environmental conditions [3].
Our results suggest that silver firs reacted to warminginduced drought in the 1980s by reducing their carbon allocation to stem growth, while needle production and retention may not have been so affected by this stressing event. However, allocation to shoot formation and perhaps root growth had likely been lowering thereafter, until the temperature rise in 1997 and 1998, which resulted in increased crown defoliation, massive decline, and enhanced mortality since 1999. Decreasing growth trends are among the most obvious growth-related characteristics of declining and dying trees [40]. The process of tree dying may take decades [13,41], and it is driven by a sequence of multiple stress factors [3,4,42,43].
Growth decline generally started several decades before tree death in the Norway spruce Swiss stands studied by [44]. Recurring drought, as judged by an objective drought-stress indicator, correlated well with modelled mortality probability for Scots pines (Pinus sylvestris) in a dry inner-alpine Predicted (logistic regression model) Figure 4: Best predictors of silver fir decline obtained by logistic regression models. Variables are plotted for declining and non-declining trees. The growth sensitivity to climate was estimated as the determination coefficient between basal area increment (BAI) and annual climatic variables (SPI, standardized precipitation index; P-ETP, annual values of the difference between precipitation-P-and potential evapotranspiration-ETP). The percentages displayed near the boxes indicate the frequency of declining (trees with percentage of crown defoliation above 50%) and non-declining trees which were correctly classified. In the box plots, error bars represent the 5th/95th percentiles, boxes indicate the standard errors, solid lines are the medians, dashed lines are the means, and points are outliers.
Swiss valley [39]. Furthermore, these authors concluded that a variable describing the most recent growth trend is needed to adequately predict tree mortality risk. The growth trends of old Pinus nigra subsp. salzmannii trees, subjected to dry climatic conditions in south-eastern Spain forests, were mostly negative over the past century [45]. Declining Slovenian silver firs showed also decreasing growth trends for about fifty years [16,46]. The growth curves of dead and living Scots pine trees at three different sites in the Swiss Valais started to diverge between ten and twenty prior to death [39]. Similar results were obtained for dying and living Abies pinsapo trees at low-elevation sites in southern Spain [13]. Silver fir growth is negatively affected by high temperature conditions of the growing season and the previous fall [47,48]. Therefore, warming-induced drought stress is likely more related to silver fir decline than changes in precipitation, since no year with extremely low precipitation was recorded in our study area since 1960. By opposite, annual temperature yielded systematically above average values since 1981. Rising temperatures increase the vapourpressure deficit and evaporation, resulting in greater water loss through transpiration [49]. Moreover, temperature rises may affect negatively carbon storage in silver fir because both the rates of carbohydrate consumption and respiration are strongly linked to temperature [50].
We conclude that the uncertainty on predicting forest vulnerability to climate warming depends on the different individual growth responses to climatic stress which may selectively alter species dynamics in sites with contrasting conditions or past successional trajectories. Our findings  Figure 5: Recent trends in basal area increment (BAI) and growth sensitivity to drought (P-ETP, difference between precipitation and evapotranspiration) for silver firs with contrasting defoliation degree; in the logistic regression declining trees correspond to the last category (trees with percentage of crown defoliation above 50%). In the box plots, error bars represent the 5th/95th percentiles; boxes represent the standard errors; solid lines represent the median; dashed lines are the mean values; and points are outliers.
emphasize the links between tree defoliation and growth sensitivity to warming-induced drought stress in water-limited populations near the southern distribution limit of Abies alba in the Spanish Pyrenees. We suggest that increasing growth sensitivity to rising air temperatures will likely portray widespread and severe, in terms of defoliation and mortality, decline episodes in drought-prone areas under a warmer climatic scenario.