Spatial variance of spring phenology in temperate deciduous forests is constrained by background climatic conditions

Leaf unfolding in temperate forests is driven by spring temperature, but little is known about the spatial variance of that temperature dependency. Here we use in situ leaf unfolding observations for eight deciduous tree species to show that the two factors that control chilling (number of cold days) and heat requirement (growing degree days at leaf unfolding, GDDreq) only explain 30% of the spatial variance of leaf unfolding. Radiation and aridity differences among sites together explain 10% of the spatial variance of leaf unfolding date, and 40% of the variation in GDDreq. Radiation intensity is positively correlated with GDDreq and aridity is negatively correlated with GDDreq spatial variance. These results suggest that leaf unfolding of temperate deciduous trees is adapted to local mean climate, including water and light availability, through altered sensitivity to spring temperature. Such adaptation of heat requirement to background climate would imply that models using constant temperature response are inherently inaccurate at local scale.

T he temporal variation in spring leaf unfolding (LU), including the effect of climatic warming, has been extensively documented [1][2][3][4][5][6] . In temperate and boreal regions, temporal variation of LU is dominated by the occurrence of warm spring temperatures 2 , typically quantified by growing degree days [7][8][9] (GDD), with LU occurring when the supplied GDD equals the required GDD. The GDD requirement of LU (hereafter GDD req ) is commonly defined as the accumulated air temperature above a temperature threshold over the preseason. GDD req varies from year to year in response to differences in the occurrence of low winter temperatures 9,10 (chilling, hereafter NCD), and constraints imposed by optimal daylength 11 . Also, the availability of water and light have been shown to affect LU and its GDD requirement [11][12][13] .
While the temporal variation in LU and its GDD req is extensively studied, the spatial heterogeneities of LU, and especially of its controls, have been much less studied. Land surface models assume that the drivers of the temporal variation of LU are also determining the spatial gradients in LU, but to our knowledge, this assumption has not been thoroughly tested. The temperature response of LU is often considered constant (spatially and temporally uniform), albeit species-dependent, in phenology models 14 . However, when applied at the regional scale, LU models were not able to accurately reproduce the observed spatial variation of LU [14][15][16][17] .
Recent warming trends have advanced LU for virtually all temperate-zone tree species, but there is large spatial heterogeneity in these advancing LU trends at the global scale 18,19 . At the regional scale, different warming trends among sites have been proposed as the main cause for these spatially differing trends in mean LU date between high and mid-latitudes 20,21 , low and high elevations [22][23][24] or coastal and inland areas 25,26 . Precipitation has also been identified as a key controlling factor of spatial differences of vegetation green-up in arid and semi-arid regions 12,27 , while at high latitudes, photoperiodic controls were proposed as an evolutionary safety mechanism that mitigates the risk of frost damage 28,29 and causes LU to occur later at higher latitudes.
Studies assessing the effects of multiple environmental variables have highlighted complex and species-specific behaviors 12,27 , mainly due to the mixing of temporal and spatial aspects of phenology, rather than trying to deconvolute these two components. Indeed, the spatial and temporal variances of phenology are expected to reflect two different mechanisms of control. First, short-term, fast responses to changing weather should drive temporal variations in LU and its GDD req , aiming to maximize fitness under inter-annually varying weather conditions. Second, an adaptive response to local biogeographical conditions may maximize tree fitness under the local long-term mean climatic conditions and would select for spatially optimized LU and climate sensitivity, inducing spatial variation therein. Biogeographic constraints on LU include all environmental variables that impose long-term adaptation of LU and its GDD req to optimize fitness under local conditions. These include climatic variables, such as site-specific occurrence of late frost events, drought occurrence, low or high light extremes, that may need to be avoided and therefore require shifts in growing season to enable maximum tree fitness. Also, site-specific interactions with neighboring competitors, pathogens and herbivores may induce spatial differences in LU and its weather dependency, in order to maximize tree fitness. Taken together, this suggests a complex response of plant phenology to climate change, but also that models of LU that apply spatially uniform parameters may not capture the observed patterns of LU and its GDD requirement.
The key hypothesis that we explore in this study is that longterm mean background biogeographical conditions determine the spatial heterogeneity of spring LU and its GDD req , reflecting evolutionary mechanisms through which plants have adjusted their growth strategies in order to maximize their fitness under those specific biogeographical conditions. We argue that biogeographic constraints on plant phenology can be detected by analyzing the spatial response of stands long-term mean LU and GDD req , instead of their inter-annual variability. In this respect, we hypothesize that not only spring, but also mean growing season conditions are important controls of the spatial differences in leaf phenology among different locations. If so, these evolutionary mechanisms control the sensitivity of LU to short-term spring temperature variations, and consequently are key components of the observed spatial variability in LU and its GDD req .
Here, we relate spatial variations in LU and GDD req in 8 dominant European deciduous tree species (see Methods) to long-term cross-sites differences in temperature, radiation and aridity (Supplementary Table 1). We combine long-term LU observations at 27790 sites over 1970-2016 (Supplementary Fig. 1) with climatic data from CRU-NCEP 30 at a spatial resolution of 0.5°( Supplementary Figs 2 and 3) to estimate the effects of long-term mean temperature, light (radiation and daylength) and aridity as determinants of the spatial variance of LU and GDD req using statistical models. Our study provides evidence for a significant control of leaf unfolding by long-term background climatic conditions across sites, potentially representing longterm adaptation of species. For all species: (1) the spatial variance of LU and GDD req supports previous correlations with temperature and chilling requirement, but is also determined by radiation intensity (W m −2 ) and site aridity, (2) the spatial variance of GDD req is better explained by preseason radiation intensity than day length at LU, (3) LU and GDD req are more sensitive to water availability on drought-prone sites and (4) LU and GDD req respond both to long-term preseason and growingseason climatic conditions. These findings suggest that at least two mechanisms influence spring phenology: (i) the direct sensing of meteorological conditions during spring to optimize the restart of plant activity and (ii) the long-term adjustment of bud sensitivity to spring meteorological conditions in order to cope with growing season pressures at sites. It also suggests that common GDD req and NCD metrics used for simulating the temporal variability of LU are not suitable for spatial studies and should be used with caution.

Results & Discussion
Importance of background climate for LU spatial variability. Among-site spatial variability of long-term mean annual LU (day of year -DoY-) and GDD req (°C; see Methods) was of the same order of magnitude than the typical interannual variation in LU and GDD req (Table 1). Figure 1 and 2 present the importance of each climatic variable (taking into account co-linearity and spatial autocorrelation; see Methods) in determining median LU date and GDD req . The full model, including long-term temperature, incident surface radiation and water availability predictors, explained 61 ± 7% of the spatial variance of LU among species (regression coefficients in Supplementary Table 2). Chilling and GDD req , which are functions of temperature and are additionally modulated by possible biotic adaptation to site-specific conditions, accounted for only half of the spatial variance of LU (Fig. 2a). Average preseason temperature (TP) was positively correlated with site LU and was selected as the best predictor instead of GDD req for half of the studied species. In addition to preseason temperature conditions, growing season temperature (TG) was found to explain 29 ± 5% of the spatial variability in LU. Considering all temperature-related variables as predictors, including preseason and growing season, models captured 52 ± 6% of total LU's spatial variance.
Long-term mean LU across sites was also significantly positively correlated with long-term growing season mean incident shortwave radiation (SWG, visible and near infrared) for all species (Fig. 2a), while no clear effect of long wave radiation (LWG, infrared) was observed. SWG contributed an additional 9 ± 2% of explained variance in LU. The importance of SWG in the model determining LU was low compared to that of temperature ( Fig. 1a), but still statistically significant (Supplementary Table 2). We also examined the relations between LU and precipitation (P), soil-moisture content (SM) and the ratio of actual to potential evapotranspiration 31 (αE; see Methods). Longterm mean LU was significantly and negatively correlated to longterm mean growing-season P (PG) for A. hippocastanum, B. pendula, F. sylvatica, F. excelsior and Q. robur, for which it captured 4.0 ± 0.6% of LU's variance. Only T. cordata exhibited a negative correlation with preseason P (PP) and SM (SMP). Growing season SM (SMG) was not correlated to LU, while αE showed a weak, but significant, positive correlation for F. sylvatica and F. excelsior, consistent with the relationship with PG ( Fig. 2a; Supplementary Table 2).
Importance of background climate for GDD req spatial variability. Long-term temperature, incident surface radiation and water availability explained 54 ± 12% of the spatial variance of GDD req among species (Figs 1b and 2b; regression coefficients can be found in Supplementary Table 3). Large differences in the impact of long-term conditions on GDD req are observed among species, with background climate explaining between 33% of the spatial variance in GDD req for B. pendula and 73% for S. aucuparia.
For all species, GDD req was negatively correlated with NCD and positively correlated with TG, together explaining 52 ± 4% of the spatial variance in GDD req (Fig. 2c). The proportion of spatial variance in GDD req explained by preseason incoming shortwave radiation (SWP) was about 30%. Day length was selected as a predictor of GDD req only for A. hippocastanum and S. aucuparia. Growing season light conditions played a minor role, albeit still statistically significant, with incoming longwave radiation (LWG) capturing about 4.0 ± 1.5% of GDD req 's explained variance (Fig. 2c). GDD req was significantly correlated with both PG and PP, capturing together around 7.5 ± 2.3% of GDD req variance. However, no significant correlations were found with soil moisture, nor αE.
Different response of phenology in drought-prone sites. While precipitation was a weaker determinant of the spatial variance in GDD req than temperature and light, for all species we observed that the spatial variance in GDD req decreased with increasing drought stress (αE) during the growing season ( Fig. 3; Supplementary Fig. 4), suggesting that frequent and/or more intense water stress during the growing-season results in adaptation of spring phenology in the long-term. Based on the observed relationships in Fig. 3, we therefore selected drought-prone sites (α E < 0.9) to assess if their responses to temperature, radiation and water differed (results for wet sites can be found in Supplementary Fig. 5 and Supplementary Tables 6 and 7).
At drought-prone sites, TP outcompeted GDD req as the main determinant of the spatial variance of LU. The contribution of temperature (NCD, GDD req , TP and TG) did not change at drought-prone sites relative to all sites (Fig. 2b). However, the relative importance of preseason precipitation (PP) was doubled at those drought-prone sites compared to the wet-sites analysis (Supplementary Tables 4 and 6). Across the drought-prone sites, LU was negatively correlated with PP for A. hippocastanum, A. glutinosa, B.pendula, F. sylvatica and Q. robur. Instead, for A. glutinosa and S. aucuparia, LU was significantly and positively correlated with PG. LU also was positively correlated with α E for F. excelsior, suggesting a differential effect of summer water availability on LU for these three species when grown at droughtprone sites. We also observed a smaller and divergent effect of incoming radiation on LU at drought-prone sites (Fig. 2b) compared to the wet-sites analysis ( Supplementary Fig. 5).
As for LU, the relative importance of temperature and incoming radiation as determinants of the spatial variance of GDD req across drought-prone sites did not differ from the other sites (Fig. 2d). However, we observed a significant contribution of water availability at drought-prone sites compared to other sites. GDD req was negatively correlated with PG and positively correlated with PP, for all species except for T. cordata. GDD req was also negatively correlated with α E for A.hippocastanum and A.glutinosa, while being positively correlated for T. cordata, which was consistent with the observed correlation with PG ( Fig. 2e, Supplementary Tables 5 and 7).
While we observed significant differences in explaining LU and GDD req spatial variance at drought-prone sites, the same approach showed no differences between warm/cold sites or low/high light sites.
GDD and NCD did not fully capture LU spatial variance. We focused on the spatial variance of spring LU and its GDD req requirement. If the determinants of the temporal variation of LU and GDD req would not fully explain their spatial variation, this would imply adaptation of spring phenology to long-term local background environmental conditions. Results revealed that temperature and the interplay between chilling and heating during winter and spring, the main determinants of the temporal variation of LU, also were important factors controlling the spatial variation of spring LU, in line with previous studies 9, 10 . However, we also showed that these chilling and heat requirement metrics only captured 30% of the total spatial variance in LU for all species (Fig. 2a, Supplementary Table 2). The positive correlation observed between LU and TP was not expected.    Higher TP implies that GDD req is reached earlier and should thus correlate negatively with LU. This positive correlation thus suggests an adaptation of LU to background temperature. Even when taking NCD, GDD req and TP into account, which all correlate with pre-season temperature, TG was still selected as the major explanatory variable of LU's spatial variance. Two hypotheses could explain this result: (1) TG is a good integrative proxy of site biogeographical constraints on LU, and potentially summarizes other variables not considered in our study, or (2) trees growing at different locations have optimized the control of their LU in response to biogeographical differences in preseason and growing season conditions 32 . In both cases, this result indicates that adaptation to long-term mean site biogeographical conditions, including growing season conditions but also a suite of biological interactions that could not be included in this study, constitutes an important evolutionary mechanism to optimize LU at that location, and must be considered in addition to commonly used preseason temperatures, which do not suffice to explain the spatial distribution of LU.
New models are needed for regional studies. In our study, longterm GDD req was strongly correlated to long-term preseason   Supplementary Fig. 6), indicating that GDD models are not independent of sites TP. With constant GDD req , higher temperature would imply earlier LU. LU indeed occurs earlier in warm areas and later in cold areas, but much less than expected based on the temporal GDD req . Trees have adjusted their GDD req to avoid too late LU in cold areas and risk damage with too early LU at warmer sites. In principle, GDD req should accounts for all warming effects, and therefore be constant unless another driver is at play or adaptation has occurred. So even if LU is partially taken as input to calculate GDD req , they would be uncorrelated. The fact that they are correlated highlights that acclimation/adaptation has occurred. The correlation between GDD req and TP and between LU and GDD req are in line with previous remote sensing studies 33,34 and strongly suggest adaptation of spring phenology to long-term site temperature. It also highlights that, even if generally useful to describe the interannual variability, GDD req in its current definition is a poor proxy of biogeographical constraints of leaf unfolding.
To emphasize this point, we additionally looked at the commonly used negative relationship between GDD req and NCD. We observed that the relationship between NCD and GDD req for LU is modulated by SWP (Fig. 4). We argue that this relationship is mainly an artifact induced by the fact that we applied uniform GDD req and NCD definition for all sites. Because trees can have different sensitivity to pre-season temperature in different regions (related to light and water availability), it implies that different GDD req and NCD definitions have to be used at the spatial level for different sites if we want to be able to simulate LU at the regional scale. Compared to current regional models using constant GDD req definition independently of the studied region, our results (Figs 2 and 4) suggested that northern sites need to have lower temperature threshold for temperature sum and/or lower GDD req thresholds than southern sites.
Since leaf unfolding determines the restart of the growing season, we expect a requirement/threshold effect of light and water conditions on LU and its temperature sensitivity. Leaf unfolding was also correlated with light and water conditions. The common choice of using GDD models (even when making GDD req dependent on NCD and day length) ignores the effects of the availability of water and light. This partly explains why GDD models typically exhibit large uncertainties when used at regional scales 14,35 . Despite their capacity to well explain the inter-annual variability of LU, the strong correlation between GDD req and long-term mean background climate observed in our study suggests that GDD models, in their current rigid parameterization, are not suitable for studying phenology at the regional scale if they do not include biogeographical constraints, especially for regions where precipitation or light are key controls of LU.
Radiation intensity matters more than day length for GDD req . Recent experimental spring phenology data have indicated that only 35% of northern hemisphere woody species relied on day length (DL) as a control of LU, and that the dependence on DL occurred predominantly at mid latitudes of the northern hemisphere 36 . Other studies have described a greater impact of daytime than nighttime temperatures on LU date 37,38 as well as differences in GDD req and NCD requirement depending on DL 11 . Day length and LU date are highly collinear ( Supplementary  Fig. 6), making it difficult to separate the impacts of day length (or preseason radiation) on LU. However, several previous studies did suggest a modulation of the temperature sensitivity of LU by light 11,[36][37][38] . In our study, growing and preseason radiation intensity received by plants were clearly identified as an important explanatory variable of GDD req 's spatial variance (Figs 1 and 2) and were more important than day length in explaining the spatial differences in GDD req .
While the effect of light intensity was clear, contributions to the spatial heterogeneity of LU and GDD req were very different between SWP and LWP. The ratio of visible to infrared light determines several plant processes, such as the state of phytochrome photoequilibrium, which controls growth rate, foliar and chloroplast development and even apical dominance 39 . Some evidence also suggests that light modulates internal hormone-regulated growth 40 and protein production in plants 41 by affecting signaling pathways of ethylene and abscisic acid, two phytohormones involved in bud set and leaf development 42 . A recent experiment found a phytochrome-mediated photoperiodic control for Fagus sylvatica 43 , while a recent review highlighted the effects of the light spectrum on spring and autumn phenology 44 . We therefore hypothesize that trees adjusted their life cycle to the average light spectrum of the site that can be directly sensed by buds and plays a direct role in enhancing or inhibiting LU, but this remains to be experimentally verified also for other species.
Light plays a key role in plant activity during the growing season by regulating photosynthesis and growth. The growing season is also the period during which buds are created. In the long-term, we expect that trees respond to growing season meteorological conditions during which the formation of buds, but also carbon reserves, can be affected, which in turn affects the sensitivity of buds to temperature in the subsequent winter and spring as already shown with Populus tremula 45 . However, differences in stand canopy openness, leaf area or even plant activity can lead to large uncertainties in the potential effect of growing season radiation on tree eco-physiology due to differences in local light regimes. As for TG, both SWG and LWG may be good proxies of background biogeographical constraints without having a direct effect on phenology.
We argue that intensity of growing and preseason incoming radiation should be included in phenological studies, not only day length as a proxy of photoperiod. More research, especially experimental studies, is needed to clearly distinguish among the effects of the light spectrum, light intensity and day length on LU and its required GDD req .
The role of aridity: The water-saving hypothesis. No theory has yet been accepted that accounts for the effect of water on spring foliar phenology in temperate forest ecosystems. Here, we propose a hypothesis on the effect of drought stress on LU dates and GDD req . LU and GDD req were more strongly correlated to water availability at drought-prone sites (Fig. 2) than at wetter sites, potentially reflecting a long-term adaptation of trees to frequent drought stress. Temperature sensitivity tends to be lower in waterlimited conditions, as indicated by the higher GDD req with increasing drought stress (Fig 2 and Supplementary Fig. 4). Organogenesis and primary growth in buds have been correlated with hydraulic architecture 46 , and previous studies have highlighted a clear effect of growing-season water stress on bud production and foliar development 47,48 . Here, we highlight an additional effect on the sensitivity of buds to temperature in spring.
We speculate that bud acclimation to previous drought may represent a water-saving strategy. By decreasing the bud sensitivity to TP, trees delay LU and the associated start of evapotranspiration 49 . The identification of preseason precipitation as an important control of the spatial variance of LU at drought-prone site is in line with this hypothesis and might also reflect a safety mechanism by which plant delay leaf unfolding until water is available for the restart of plant activity. As a longterm strategy, a delay in evapotranspiration will lead to a slower depletion of water resources at the beginning of the growing season and reduces the risk of water stress during late spring and summer when radiation is more favorable for photosynthesis.
Species show different responses to long-term constraints. All species exhibited the same response to temperature, however a few species responded differently to water availability at drought-prone sites. A. glutinosa, T.cordata and F.excelsior showed opposite correlation with water availability compared to the other species (Fig. 2). T. cordata and F. excelsior were proposed as species with a relatively high drought tolerance compared to other European species 50,51 , while A. glutinosa naturally occurs in wet sites. This different behavior observed for A. glutinosa, T.cordata and F. excelsior could reflect either their drought tolerance, or the fact that differences in edaphic heterogeneity were not captured by aggregating data at the pixel level, resulting in a mismatch between actual water availability and the estimates of αE used here.
Differences were also observed for some species regarding their heat requirement. Biogeographical conditions of temperature, light and water only captured 33, 44 and 49% of GDD req spatial variance for B.pendula, F. sylvatica and A. glutinosa, respectively. Previous studies showed a reduced sensitivity of leaf unfolding to climate warming in the last decades, mainly attributed to plant plasticity 52 . However, how trees acclimate or adapt to future climate change remains unclear and might be species-dependent, due to differences in temperature, light and water sensitivities.
The remaining unexplained variance of LU and GDD req can be attributed to uncertainties in observations and climate data, with a potential effect of data aggregation at the pixel level, but also of unaccounted biotic and edaphic factors (e.g., stand age, biogeography of pathogens or mycorrhizal associations, soil structure and fertility, etc.).
In the end, the co-limitation of spring phenology by light and aridity may account for why LU does not keep pace with climate change 19,52 , which may have vast and far implications on the carbon cycle 53 , with a possible alteration of the competitive balance among species 54 .
Tree seasonality affects spring phenology. Temperature, availability of light and water and their interactions thus affect the spatial heterogeneity in LU date and its associated GDD req . On one hand, our results showed that LU date adjusted to spring meteorological conditions, for which trees are co-limited by TP, SWP, LWP and PP. On the other hand, we observed different responses of spring phenology to preseason and growing-season meteorological conditions, especially visible for drought-prone sites, highlighting that tree adjust their phenology to cope with seasonality on the long-term.
Since buds are formed during the growing season, we argue that any effect of growing season pressure might in turn have a feedback on tree seasonality, with a potential long-term adjustment of ecophysiological responses to these constraints. In the long term, it may also alter the restart of plant activity in spring and thereby optimize the long-term growth, reproduction, or survival of the trees, which influence the restart date of physiological processes.
The impact of elevated air temperature during the growing season has been proposed to affect spring leaf unfolding by modifying growth cessation and dormancy induction in temperate and boreal trees 55 . Drought stress also affects the onset of senescence and possibly dormancy 56,57 . Precipitation has indeed been shown to play a large role in vegetation dynamics during the senescence period of deciduous forests in the Northern Hemisphere 58 . A delayed dormancy induction can translate into a delayed leaf unfolding during the next spring via the later start date of chilling and GDD accumulation 59 . Our results are fully in line with these observations and provide additional evidence that the timing of phenological events is impacted by previous phases in the annual cycle of trees.
We argue that tree seasonality and long-term biogeographical constraints are too often overlooked and should be taken into account in phenological studies. Further research on the effect of tree seasonality on inter-annual variability of phenology is needed to clearly identify the role played by biogeographical constraints.
In conclusion, assessing the long-term spatial variance of LU and GDD req is a step in developing a unified framework that will allow an understanding of the multiple control of climate on plant phenology. Future research on the importance of plant phenology on ecosystem functioning should focus on space-time interactions with environmental conditions specifically to address: 1) the effects of light and aridity on bud sensitivity to temperature, and 2) the potential coordination between plant processes and phenology that could account for a co-limitation by temperature and the availability of light and water. In line with these recommendations, the use of current, constant, GDD req and NCD metrics for the study of spatio-temporal patterns in plant phenology should be used with caution.
Daily temperature, precipitation and incoming radiation (both shortwave [visible and near infrared] and longwave [infrared]) were retrieved from the CRU-NCEP climatic data set 30 at a spatial resolution of 0.5°. Day length was calculated for each site using the 'geosphere' R package 60 . Monthly SM (10-40 cm) was retrieved from the land data assimilation systems (GLDAS) data set (https://ldas. gsfc.nasa.gov/gldas/) aggregated at a spatial resolution of 0.5°to be consistent with the CRU-NCEP data.
We approximated drought conditions using the ratio of actual to potential evapotranspiration as a proxy for drought periods (αE). This ratio accurately represents droughts 61 and was calculated using the SPLASH model 31 . αE was calculated daily, but only the growing-season αE was used to calculate long-term drought stress.
Analyses. We corrected site temperature for altitudinal differences between the site and the mean elevation of each grid cell of the CRU-NCEP data set 14 using a gradient of 6.4°C km −1 and then estimated the heat requirement and associated chilling for each LU observation. Heat requirement (GDD req ) corresponds to the sum of mean daily temperatures above a threshold of 5°C calculated from 1 January to the date of LU and was calculated for each observation as: GDD req defines the heat requirement of buds at the observed LU date, estimated as the accumulated daily air temperature (GDD d ), where T d is the mean daily air temperature, T th the temperature threshold for GDD accumulation (5°C) and t o the starting date (1 January).The number of chilling days (NCD) was calculated as the number of days from 1 November to the LU date with mean daily temperatures between 0 and 5°C.
We were interested in the spatial variance of LU and the heat requirement, so we used the median LU date and the corresponding median GDD for 1970-2016 for each site, assuming that the medians represented the long-term "optimal" LU date and GDD for the background climatic and soil conditions. We only retained sites with more than 5 years of observations and removed years with a LU date outside two interquartiles around the median distribution (i.e., outside days 80-152), which potentially represent a response to extreme events. The same analysis was performed using sites with more than 10 years of observations and led to similar results. We analyzed the correlations between biogeographical variables, LU dates and GDD for eight species of dominant European deciduous trees with many records in the PEP725 database: Aesculus hippocastanum, Alnus glutinosa, Betula pendula, Fagus sylvatica, Fraxinus excelsior, Quercus robur, Sorbus aucuparia and Tilia cordata.
Long-term climatic variables (i.e., over 1970-2016) were averaged for two periods of the year: the growing season, defined as the period between days 180 and 250, and the preseason, defined as the three months before LU. Since we are interested in biogeographical constraints, not the temporality of processes, and because no information about the length of the growing season was available in the PEP data set, selecting a constant summer period ensures that we have a representative period for all trees and all years. It allowed a consistent comparison between sites without introducing bias induced by different growing season lengths. This period also corresponds to the peak of plant activity in temperate ecosystems, and over which we are more likely able to gather information about water and temperature pressures of each site.
As for LU, extreme climatic years were excluded as we seek to estimate the average response of the vegetation. We then analyzed the spatial relationships between long-term climatic variables (averaged over the different periods of the year), LU and GDD. We proceeded in four steps: 1. We first assessed potential collinearity issues between variables by examining pairwise correlation coefficients 62 (Supplementary Fig. 6); 2. We then selected relevant predictors of observed LU and GDD req using penalized elastic net regressions (glmnet function from the glmnet 63 package) in combination with collinearity information from step 1). Interaction terms were not included in the analysis. In step 2-4) predictors were all standardized in order to represent the relative contribution of each variable in explaining LU and GDD variability ( Supplementary Tables 8 and 9).
3. After selecting relevant variables in step 2, we assessed the remaining collinearity between variables by estimating their respective variance inflation factor (VIF) as: where the VIF for variable j is the reciprocal of the inverse of R² from the regression. VIF values increase with collinearity, and arbitrary threshold of 5-10 are commonly used to define high VIF values. Here, when two predictors exhibited potential collinearity (i.e., high VIF values), we removed the one with a VIF value higher than 4 and the lowest correlation coefficient with LU or GDD using the stepVIF function from the VIF 64 package. 4. We finally assessed the spatial structure of residuals of the reduced model from step 3) with semi-variograms ( Supplementary Figs 7-9). We performed generalized least square regressions taking into account the spatial structure of residuals to correct coefficients using the gls function from the nlme 65 package. Different spatial structures were tested (linear, exponential, spherical, Gaussian and rational quadratic) and the best model was selected using AIC criterion (Supplementary Table 10).
All the above steps were applied to each species separately. Because results were consistent between species (Fig. 2) regarding aggregation uncertainties and because we aimed at exploring biogeographical constraints at the regional scale for all deciduous tree species, we applied the same approach to the full dataset with all species pooled together. Analyses were performed with the R v.3.5 software 66 .
Reporting summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.