Next Article in Journal
Developing a Method to Extract Building 3D Information from GF-7 Data
Next Article in Special Issue
Exploring the Applicability and Scaling Effects of Satellite-Observed Spring and Autumn Phenology in Complex Terrain Regions Using Four Different Spatial Resolution Products
Previous Article in Journal
Evidence of Climate Change Based on Lake Surface Temperature Trends in South Central Chile
Previous Article in Special Issue
Quantification of Urban Heat Island-Induced Contribution to Advance in Spring Phenology: A Case Study in Hangzhou, China
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Remote Sensing
Volume 13
Issue 22
10.3390/rs13224538
remotesensing-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Specific Drivers and Responses to Land Surface Phenology of Different Vegetation Types in the Qinling Mountains, Central China

by
Jiaqi Guo
1,
Xiaohong Liu
1,*,
Wensen Ge
1,
Xiaofeng Ni
2,
Wenyuan Ma
3,
Qiangqiang Lu
1 and
Xiaoyu Xing
4
1
School of Geography and Tourism, Shaanxi Normal University, Xi’an 710119, China
2
Key Laboratory for Earth Surface Processes of the Ministry of Education, Institute of Ecology, College of Urban and Environmental Sciences, Peking University, Beijing 100871, China
3
School of Science, University of New South Wales, Canberra, ACT 2600, Australia
4
Qinling National Botanical Garden, Xi’an 710061, China
*
Author to whom correspondence should be addressed.
Remote Sens. 2021, 13(22), 4538; https://doi.org/10.3390/rs13224538
Submission received: 9 October 2021 / Revised: 4 November 2021 / Accepted: 8 November 2021 / Published: 11 November 2021
(This article belongs to the Special Issue Remote Sensing of Land Surface Phenology)
Browse Figures
Versions Notes

Abstract

:
Land surface phenology (LSP), as a precise bio-indicator that responds to climate change, has received much attention in fields concerned with climate change and ecology. Yet, the dynamics of LSP changes in the Qinling Mountains (QMs)—A transition zone between warm-temperate and north subtropical climates with complex vegetation structure—under significant climatic environmental evolution are unclear. Here, we analyzed the spatiotemporal dynamics of LSP for different vegetation types in the QMs from 2001 to 2019 and quantified the degree of influence of meteorological factors (temperature, precipitation, and shortwave radiation), and soil (temperature and moisture), and biological factors (maximum of NDVI and middle date during the growing season) on LSP changes using random forest models. The results show that there is an advanced trend (0.15 days/year) for the start of the growing season (SOS), a delayed trend (0.24 days/year) for the end of the growing season (EOS), and an overall extended trend (0.39 days/year) for the length of the growing season (LOS) in the QMs over the past two decades. Advanced SOS and delayed EOS were the dominant patterns leading to a lengthened vegetation growing season, followed by a joint delay of SOS and EOS, and the latter was particularly common in shrub and evergreen broadleaved forests. The growth season length increased significantly in western QMs. Furthermore, we confirmed that meteorological factors are the main factors affecting the interannual variations in SOS and EOS, especially the meteorological factor of preseason mean shortwave radiation (SWP). The grass and crop are most influenced by SWP. The soil condition has, overall, a minor influence the regional LSP. This study highlighted the specificity of different vegetation growth in the QMs under warming, which should be considered in the accurate prediction of vegetation growth in the future.

Graphical Abstract

1. Introduction

Vegetation phenology is the seasonal timing of lifecycle events, such as leaf emergence, flowering, leaf coloration and fall, and it has become an important topic in the field of climate and ecology as a sensitive and precise indicator that is responsive to climate warming [1,2]. Shifts in spring and autumn vegetation phenology caused by climate warming can differentially alter the length of the growing season, which affects carbon, water, and energy exchange between terrestrial ecosystems and the atmosphere [3,4,5]. Recent studies have reported that in addition to climatic factors, soil and biological factors also influence shifts in vegetation phenology by affecting plant growth processes in the context of ongoing global climate change [6,7], due to the poor interpretation of phenology shifts among different vegetation types [8,9]. Hence, it is essential to study the dynamics and drivers of phenology among different vegetation types to improve phenology models and enrich our understanding of the carbon cycle of terrestrial ecosystem.
With the application of remote sensing in monitoring vegetation phenology, we traditionally use the term land surface phenology (LSP) to denote the dynamic variations in vegetation land surface as observed from satellite imagery [10]. Satellite-derived LSP metrics are usually focused on the start (SOS) and end (EOS) of growing seasons [11]. Satellite-based studies have shown that SOS was advanced by 10.6 days (i.e., 5.4 days per decade) throughout Europe and by 14 days (i.e., 7.9 days per decade) in temperate China before 2000 [12,13]. However, this trend of SOS advancement may have slowed or even reversed since the 2000s. For instance, in the entire northern hemisphere, it was advanced by only 0.2 days during 2000 to 2008, but a delayed SOS was revealed in the Tibetan Plateau [14,15]. Regarding satellite-derived EOS, the published results have not always been consistent. Across the entire Northern Hemisphere, EOS was delayed at a rate of 2.2 days per decade during 2000–2008 [14]. In the Yellow River Basin, EOS was delayed by 5.6 and 3.4 days in 1982–1999 and 2000–2015, respectively [16]. In the Qinghai-Tibet Plateau, however, Wang et al. [17] reported the opposite phenology change trends, in both the east and west zones. Overall, these varied results might be due to different study areas, periods, and methods of extracting phenology metrics. However, few studies have focused on the diversity of phenology across different vegetation types. In particular, the dynamic phenology characteristics of herbaceous or shrubs and evergreen forests in the Qinling Mountains have been little studied.
To date, the processes and drivers governing LSP remain poorly understood. Several studies reported that temperature is a major driver of early spring leaf development and delayed autumn leaf fall in plants and has less control over autumn phenology than spring phenology [18,19]. The impact of precipitation on LSP processes is mainly directed at plants in arid and semiarid regions, where water deficits limit the use of light and heat conditions by plants in arid and semiarid areas [20]. Some studies further considered solar radiation (i.e., shortwave radiation) and found that increasing photosynthetic active radiation can promote earlier leaf germination and delaying leaf senescence [21]. Besides meteorological factors, soil factors and biological factors have also been shown to be drivers of LSP change processes [22,23,24]. Soil temperature and moisture information, due to the prevalent freezing—Thawing process of soil in alpine and arctic regions, could more directly control vegetation growth; for example, soil wetting will, to some extent, reduce the effect of soil warming on LSP changes [22,25]. Slight fluctuations in the time interval between the middle date (MD) of the growing season and the autumn phenology have a strong effect on regulating EOS [7]. Peak growth in summer (i.e., maximum NDVI during the growing season, MN) can have an impact on vegetation greening and senescence, and its unique vegetation growth patterns may result in different allocations of green or carbon across the growing season [24]. To date, how these meteorological, soil, and biological factors affect regional-scale LSP variations has not been clearly and consistently studied, which seriously affects our ability to predict LSP periods.
The responses of vegetation growth processes to climate change are inherently nonlinear [26]. Random forest (RF), as a nonparametric multivariate method, can explain nonlinear processes to a large extent [27]. The advantage of the RF model is that it can consider many predictor variables and nonlinearly determine the relative importance of each predictor variable [28]. Due to its high efficiency in handling the potentially complex relationships between LSP periods and meteorological, soil, and biological factors, RF has been widely and successful applied in recent ecology studies [6,29].
The Qinling Mountains (QMs), rich in vegetation types, represent a demarcation line of climate in China and are also an area characterized by sensitivity or response to climate change, with a significant upward temperature trend in the past half century [30]. Here, we used the satellite derived normalized difference vegetation index (NDVI) records (2001–2019) from MOD13A2 to extract the LSP dates of QMs. The objectives were (a) to explore the temporal and spatial trends of LSP in the QMs, (b) to quantify the relative contribution of SOS and EOS of different vegetation types to the length of growing season (LOS) and to determine the dominant growth pattern during the growing season, and (c) to simulate the LSP dates and assess the relative importance of meteorological, soil, and biological factors on the interannual variations in LSP. This study focuses on the specificity of different vegetation growths in the QMs, and the results are helpful for future accurate prediction of vegetation growth and to develop scientific management strategies.

2. Materials and Methods

2.1. Study Area

Qinling is the highest mountain range in the central region of China and also a geographical boundary between the north and the south of China, with an elevation range of 51 to 5120 m and a spatial range from 30.8° to 35.5°N to 102.5° to 114.6°E (Figure 1a). Its climate differs significantly between north and south, with a humid northern subtropical climate in the south and a warm temperate semi-humid and semiarid climate in the north (Figure 1a). It has also been classified as one of the critical terrestrial biodiversity areas of world significance [31], with mixed coniferous and deciduous broadleaved forests widely distributed on its northern slopes, while the southern slopes are dominated by mixed evergreen and deciduous broadleaved forests (Figure 1b).

2.2. Datasets

NDVI is the most commonly applied vegetation index to characterize vegetation greenness and is strongly correlated with vegetation photosynthetic activity [32]. Climate change is the main factor affecting the change of vegetation greenness, and this change can be reflected by the spectral information of NDVI images. In this study, we used the NDVI datasets generated from NOAA/AVHRR series satellite images by the NASA MODIS13A2 group (Table 1). We used this dataset to extract the LSP dates for QMs from 2001 to 2019. We also excluded areas of bare soil/sparse vegetation with an annual average NDVI of less than 0.1 [33].
Meteorological data, including daily mean air temperature, daily total precipitation, and daily mean shortwave radiation and soil data, including daily mean soil temperature and moisture in 0–100 cm soil layer, from 2001 to 2019 were used in this study. These gridded data were derived from the ERA5-Land hourly data (Table 1). Moreover, we transformed the hourly climate data (24 hourly data were averaged for temperature, solar radiation, soil temperature and humidity, 24 hourly data were summed for precipitation) to daily-scale temporal resolution and resampled meteorological data to the same resolution as MODIS13A2 data. A time lag of 30 days before SOS and 60 days before EOS is defined as preseason.
The 300 m spatial resolution Climate Change Initiative Land Cover (CCI-LC) maps from 2001 to 2019 were available from the European Space Agency (ESA) (Table 1). CCI-LC discriminates 22 classes of land cover. In this study, we resampled these maps to 1 km and analyzed only pixels of unchanged vegetation types containing evergreen needle leaved forest (ENF), evergreen broadleaved forest (EBF), deciduous broadleaved forest (DBF), mixed forest (MF), shrubland (SL), grassland (GL), and cropland (CL).

2.3. Retrieval of Phenology Metrics from NDVI Time Series Data

The premise of quantitatively analyzing the phenology changes is to derive several key phenology metrics: SOS, EOS, LOS, MN, and MD (Figure 2). In this study, we firstly stacked the NDVI images from 2001 to 2019 in chronological order and smoothed the NDVI time series with a Savitzky-Golay (SG) filter for each pixel per year. The SG filter was chosen because it can best preserve the temporal vegetation dynamics and minimize atmospheric contamination and has also been integrated into the processing of the MODIS phenology product [34]. The smoothed data was used further for extracting phenology metrics of different vegetation types by detecting the inflection point (i.e., date) when the NDVI time series begins to ascend or descend for the specific year. This is the derivative method which the phenology metrics were extracted for each pixel per year, whereby the maximum value of NDVIratio corresponds to the greatest change of the smoothed NDVI time series [35]. Equation (1) is given as
NDVI r a t i o ( t ) = NDVI t + 1 − NDVI t NDVI t
where NDVI r a t i o t is the calculated relative changing rate of NDVI at time t and NDVI t is the NDVI value at time t. Occurrence dates were obtained using these smoothed NDVI time series. SOS and EOS dates were determined as the day with the maximum and minimum NDVI r a t i o . LOS was determined to be the difference between EOS and SOS. MN was defined as the peak of vegetation growth, i.e., the NDVI value corresponding to NDVI r a t i o closest to zero. MD date was the middle date between the EOS date and the SOS date. The description of phenology metrics correlations is shown in Figure 2.

2.4. Method and Statistical Analysis

2.4.1. Trend Analysis

The method used in this study is shown in Figure 3. The spatiotemporal trends of LSP during 2001–2019 were estimated using Sen’s slope method, also known as the Theil-Sen median method. The method is a robust nonparametric statistical method for trend calculation that is insensitive to measurement errors and is far more accurate than nonrobust simple linear regression [36]. Sen’s slope was calculated using Equation (2):
S e n ’ s   s l o p e = M e d i a n x j − x i j − i
where the median is the mean value of all the slopes, and xi and xj represent the LSP dates of years i and j. A negative Sen’s slope indicates an advancing trend, whereas a positive Sen’s slope indicates a delaying trend.
Then, we used the Mann-Kendall (MK) method to test the significance of time series trends, which is a nonparametric statistical test and is robust to outliers [37]. We used the normal cumulative distribution function to determine the p-value of the MK test statistic with a significant confidence level of p < 0.1. In this study, we used Sen’s slope and MK test to trend analysis and significance test the spatial distribution (average growth) and interannual variation (interannual growth) of LSP for different vegetation types from 2001 to 2019, all using MATLAB 2017a were completed.

2.4.2. Change Pattern and Relative Attribution Analysis

To further understand the seasonal changes of vegetation growth in the QMs in the past two decades, we divided the trend changes of LSP (SOS, EOS, and LOS) into six combinations, where each combination represents a pattern of vegetation growth (Table 2). We used spatial analysis to count the proportion and significance of different vegetation for each pattern and to derive the dominant pattern of seasonal changes in the growth of each vegetation. All analysis was accomplished using ArcGIS 10.4.
To evaluate the symmetry of SOS and EOS for LOS changes, we used the C-index proposed by Garonna et al. [38] to calculate the relative contribution of trends in SOS and EOS for LOS changes. It was calculated as follows:
C = − 1 + 2 ⋅ S O S s l o p e S O S s l o p e + E O S s l o p e
where SOSslope and EOSslope are the Sen’s slope of SOS and EOS, respectively. A positive C value indicates that the trend in LOS is mainly attributable to changes in EOS, and a negative C value indicates that the trend in LOS is mainly attributable to changes in SOS. The variation of C value is from −1 to 1.

2.4.3. Analysis of the Relative Importance of Different Drivers

Based on previous studies on the drivers of interannual variation in vegetation phenology [18,19,20,21,22,23,24], we selected drivers such as Table 3 to simulate SOS and EOS. These drivers are divided into three main categories: meteorological, soil and biological factors, which are further divided into preseason cumulative values and cumulative values throughout the growing season.
In this study, the RF model was used to assess the relative importance of the drivers affecting interannual variations in LSP. First, for each pixel, we calculated the partial correlation coefficients between environmental factors (TP, TG, PP, PG, SWP, SWG, STP, STG, SMP, SMG, MN, and MD) and SOS and EOS during 0, 1, 2, … n months before SOS and EOS, and we separately derived the highest correlation with SOS and EOS for the time range of 30 days before SOS and 60 days before EOS for environmental data. Moreover, a subset of variables highly correlated with SOS and EOS was selected, and the values of the variables at selected years and locations (spatiotemporal models) were combined into a set of input feature vectors that are used as inputs for the RF algorithm. Then, these feature vectors were divided equally into two subsets, with 2/3 of the dataset used for model training (in bag) and the remaining 1/3 of the dataset used as an additional test of the RF internal computation (out of bag, i.e., OOB) to estimate the importance of each variable. Variable importance can also be measured by OOB, which compares the increases in OOB error with that variable randomly permuted and all others unchanged [39]. The importance score of a variable is as follows:
V I X j = 1 n t r e e ∑ t e r r ’ O O B t j − e r r O O B t j
where X j is the jth variable, ntree is the number of trees, e r r O O B t j is the OOB error of each tree t, and e r r ’ O O B t j is the OOB error when X j is permuted, while all other variables remain unchanged among OOB data. For regression, the OOB error is the mean square error.
Finally, to optimize the model, the hyperparameter search was used to select the best tested hyperparameter set, and the optimal model was trained on the whole training set for accurate prediction of the LSP dates. Besides modeling multi-year-scale LSP variation for the entire region, subregional variation according to different vegetation types was also modeled to quantify the relative importance of different drivers. These models were constructed using the RF package in R statistical software.
To evaluate the predictiveness of the model and to further test the applicability of the model, we used a randomly selected subset (1/3 of the dataset) for model validation. Both the proportion of explained variance (R2) and the root mean square error (RMSE) were used to assess the performance of the model on the complete datasets. The stability of the model fit is explained by the mean absolute error (MAE) [40]. These statistics were calculated as follows:
R 2 = 1 − ∑ i = 1 n ( y i − y ⌢ i ) 2 ∑ i = 1 n ( y i − y ¯ i ) 2
R M S E = ∑ i = 1 n ( y i − y ⌢ i ) 2 n − 1
M A E = ∑ i = 1 n y i − y ⌢ i n
where y i is the observed satellite-based SOS or EOS values, is y ⌢ i the predicted RF-based SOS or EOS values, and y ¯ i is the mean satellite-observed SOS or EOS values of all selected test pixels for 2001–2019. n is the sum of all selected test pixels.

3. Results

3.1. Spatiotemporal Variations of Phenology Metrics in the QMs

Figure 4 shows the spatial variation of the annual mean LSP and their corresponding standard deviations (Std) over the study period 2001–2019. Earlier (<90 days) sites of SOS (14.4%) were located at low elevations in the central QMs, and later (>130 days) sites (6.8%) were located at high elevations in the western QMs (Figure 4a). The earliest occurrence of SOS was for SL, with a mean SOS of 97 ± 14 days, and the latest occurrence was for GL, with a mean SOS of 109 ± 12 days (Figure 4a). The Std of SOS has significant spatial variation, with larger areas located in the southwestern QMs (57.2%) having higher Std (>9 days) (Figure 4b). The overall spatial variation in multi-year average EOS was not significant, and in the southeastern QMs, EOS was mainly concentrated in 290–300 days, accounting for 43.0% of the entire study area (Figure 4c). The earliest EOS occurred in GL, with a mean EOS of 287 ± 11 days, and the latest occurred in EBF, with a mean EOS of 295 ± 12 days (Figure 4c). There was higher Std (>13 days) in the southern QMs (25.2%) compared with the northern QMs (Figure 4d). LOS had clear spatial differences, with the central QMs (15.5%) having the longest LOS (>210 days) and the western high-altitude areas (8.5%) had the shortest LOS (<150 days). For SL in the southern QMs (27.9%), its mean LOS was 193 ± 20 days and the std (>18 days) was also the largest (Figure 4e,f).
We also characterized the spatial distribution of LSP trends for different vegetation types from 2001 to 2019 (Figure 5). For the whole QMs, SOS was advanced in 67.8%, the average rate of advance was 1.5 days/decade, and 27.5% of the area (mostly located in the northern QMs) was significant (Figure 5a,b). DBF advanced at a rate of 1.9 days/decade and was the fastest compared to other vegetation types (Figure 5a). EOS was delayed in 72.1% of the region and significant for 42.1% of the region (mostly located in the southern QMs), with an average delay rate of 2.4 days/decade across the region (Figure 5c,d). EBF had the fastest delay rate of 3.3 days/decade (Figure 5c). The average rate of LOS lengthening across the study area was 3.9 days/decade, and 74.6% of the areas (mostly in the southwestern QMs) were lengthened (Figure 5e). The rate of LOS lengthening was 4.7 days/decade for EBF, fastest among the seven vegetation types. Of these areas that were lengthened, the change was found to be significant in 40.3% (mostly in the western QMs) of cases (Figure 5f).
Data of interannual variation trends and the significance of LSP for different vegetation types are shown in Figure 6. Overall, the Sen’s slope of SOS is −0.09 days/year from 2001 to 2019, but this advance is insignificant. There was a trend of significant SOS advancement for DBF and GL, at 0.16 and 0.13 days/year, respectively. EOS shows a significant delay trend over the entire region of 0.29 days/year. The trend of EOS delay was more significant for both EBF and CL compared to other vegetation (p < 0.05), and the rate of EBF delay was the fastest (0.37 days/year). LOS is significantly lengthened at a rate of 0.48 days/year. ENF, EBF, DBF, MF, and CL show a more significant trend for lengthened LOS (p < 0.05), with the fastest being for EBF (0.68 days/year).

3.2. Change Pattern of LSP and Relative Attribution Analysis

For all vegetation types, the dominant change pattern of the growing season was Type I, with a proportion of 48.4%, which implies that most vegetation in the QMs had an advanced SOS, delayed EOS, and lengthened LOS (Table 4). Type V showed the second largest proportion (15.2%) which meant that there were also many plant species on the QMs having delayed SOS, delayed EOS, and lengthened LOS. Types III, IV, and VI had the smallest proportions (all lower than 10.0%), which indicated that the probabilities of shortened LOS were very low for all plants on the QMs. For five of these vegetation types (ENF, EBF, MF, SL, and CL), the dominant change pattern of the growing season was Type I, followed by Type V. This indicates that these types of plants on the QMs had delayed EOS and lengthened LOS. The main change pattern for DBF and GL was also Type I, with Type II being the second most prevalent. This implies that there is some DBF and GL showing advanced SOS, advanced EOS, and lengthened LOS.
Figure 7 shows the significance of each pattern. For all vegetation types, lengthened LOS, advanced SOS, and delayed EOS were significant in terms of Type I, II, and V change patterns, respectively. Types I, II, and V were the top three patterns in terms of percentage, as seen in Table 4, which were also the three patterns of LOS lengthening. In the Type I, lengthened LOS is significant for most vegetation types (EBF, MF, SL, GL, and CL). In terms of Type II change, the SOS of ENF, DBF, MF, GL, and CL were all significantly advanced, and advanced SOS resulted in lengthened LOS. Fewer changes in LSP trends were significant in terms of Type III, IV, and VI changes, only SL and GL showed significant Type IV changes, and SL (15.5%) and GL (12.2%) accounted for a large proportion of Type IV changes. In terms of Type V, the EOS of all six vegetation types (ENF, EBF, DBF, MF, GL, and CL) was significantly delayed. Delayed EOS resulted in lengthened LOS, a pattern which also happens to be the second largest in terms of proportion, and this pattern is also the one we should be concerned about.
For the entire study area, the calculated C-index values were negative for 106,385 pixels and positive for 112,514 pixels (Table 5 and Figure S1). Pixels with C-index values less than 0 are mainly distributed in the easternmost and southernmost parts of the QMs, which indicates that their LOS variations are mainly controlled by SOS shifts. All other regions have pixels with values over 0, which indicates that they are controlled by EOS shifts (Figure S1, in Supplementary Materials). The percentage of LOS changes controlled by SOS and EOS is 48.6% and 51.4%, respectively (Table 5). This also shows that LOS trends, except for DBF and GL, were mainly controlled by the shift in EOS for each vegetation type. The percentage of LOS changes controlled by SOS shifts was 53.4% and 52.0% for DBF and GL, respectively. The largest percentage of changes in LOS of all vegetation types controlled by EOS was 58.3% (SL), and the smallest was 46.6% (DBF).

3.3. Drivers of Interannual Variations in LSP

For the QMs, different drivers affect the interannual variability in SOS and EOS (Figure 8 and Table 6). The SWP, MD, and STP are the three most important factors influencing the interannual SOS variation, and the relative importance accounts for 54.4% of total (Figure 8a and Table 6). The total percentage of TG, PP, TP, STG, and PG was 39.2%, and the effect of TG and PP on the interannual SOS variation was almost the same. The remaining four variables (SWG, SMG, SMP, and MN) have a very small effect on interannual SOS variation, and their combined percentage was only 6.4%. Figure 8b and Table 6 also show that SWP, PP, and MD are the three most important factors influencing the interannual EOS variation, with a total relative importance of 54.0%, and the influence of SWP is much stronger than that of PP and MD. The effects of TP, SWG, STP, PG, TG, and STG on the interannual EOS variation totaled 41.9%, and the effects of TP and SWG on the interannual EOS variation were not very different, with a relative importance of 10.2% and 9.9%, respectively. There was also little difference in the relative importance of PG, TG, and STG. The effect of the single variable of STP on the interannual EOS variation (7.2%) is much larger than the sum of SMP, SMG, and MN (4.1%).
The drivers influencing interannual variations in LSP of different vegetation types were assessed (Figure 9 and Table 6). Figure 9a and Table 6 show that the main drivers affecting the interannual SOS variation of ENF, MF, and GL were SWP, MD, and STP, and the relative importance of SWP in these three vegetation types was ranked as GL (28.6%) > MF (22.1%) > ENF (21.2%). The interannual SOS variation of EBF, DBF, and CL was mainly influenced by MD and SWP. The main drivers influencing the interannual SOS variation of SL were STP, TP, and SWP, and STP (29.2%) was the most important factor influencing the interannual SOS variation of SL. As shown in Figure 9b and Table 6, in terms of the interannual EOS variation, SWP had the strongest effect on GL (36.4%) and the slightest effect on SL (20.5%). The main drivers of interannual EOS variation of ENF, DBF, and CL are SWP, TP, and PP. SWP, PP and MD are the main drivers of interannual variations in EOS for EBF, MF and SL. Besides SWP, which is the most important driver, PP is the second factor affecting MF and SL, and MD is the second factor affecting EBF with a relative importance of 14.5%. Moreover, the main factors affecting the interannual EOS variation of GL are SWP, SWG, and TP, and the relative importance of SWG and TP is 11.5% and 11.2%, respectively.
We used 19 years of data to construct an RF model for seven vegetation types (Table S2), and we randomly sampled 1/3 of the pixel dataset to assess their linear relationships (Figure S2). The actual and predicted SOS displayed good linear relationships, with their correlation coefficient R2 values ranging from 0.900 to 0.938, RMSE values ranging from 4.90 to 6.52, and MAE values ranging from 3.93 to 5.15 (Figure S2a). Both the actual and predicted EOS also show good linear relationships, with their correlation coefficients R2 values ranging from 0.911 to 0.942, RMSE values ranging from 4.23 to 5.63, and MAE values ranging from 3.26 to 3.74 (Figure S2b). These results indicate that it is appropriate to use RF models to analyze interannual variations in LSP in the QMs.

4. Discussion

4.1. Dynamics Changes in LSP in the QMs

Understanding the interannual variations in vegetation phenology and its trends is important for recognizing the patterns of vegetation growth dynamics as a response to climate warming. Our study showed that there is an advanced trend (1.5 days/decade) for SOS, a delayed trend (2.4 days/decade) for EOS, and an overall extended trend (3.9 days/decade) for LOS in the QMs during 2001–2019. In comparison with previous studies on phenology changes, there were different degrees of advanced SOS, delayed EOS, and lengthened LOS in different study areas and periods [32,41]. For example, during 1982–2006, the SOS was advanced by 0.56 days/decade, while the EOS delayed trend rate was 5.5 days/decade, and the growing season was significantly longer by 6.06 days/decade in North America [42]. In the Tibetan Plateau region, SOS was advanced at a rate of 0.17 days/decade, EOS was delayed at a rate of 5.29 days/decade, and LOS was lengthened at a rate of 5.46 days/decade for the period 1981–2017 [41]. These results are not entirely consistent with those of studies conducted for the QMs, and the differences may be mainly due to their different target periods and the different methods of phenology extraction. However, other investigations reported that, compared to 1982–1999, the phenology trend slowed down in 2000–2008 and the changes were not highly significant [14,43]. Meanwhile, our results also show a slower change in phenology trends over the last 20 years in the QMs, and the magnitude of SOS advance is also smaller than that of EOS delay. This observation is similar to the results of Wang et al. and Xia et al. [33,44], which indicates that the satellite-observed phenology change rates slowed down during a global warming hiatus between 1998 and 2012.
We also found that the trend in phenology changes of different vegetation is also highly variable, which is related to the microclimate of different regions or the geographical variation of plant origin [45]. Our results show a significantly greater trend of SOS changes in SL than other vegetation types, indicating that the earlier the SOS, the more significant the trend in SOS variations, and that this difference may be related to plant pollination type, life type, phylogenetic and wood type, etc. [46]. However, the trend of SOS changes in CL was weaker than in other vegetation types, and this trend was insignificant because farmers controlled the sowing time in each year, resulting in significantly smaller variability in crop phenology than in field observed plants [47]. The results show that the trend of delayed EOS is significantly stronger for EBF than other vegetation types, mainly because EBF is mainly located in the region south of the QMs, with a humid northern subtropical climate. Some researchers have shown that in subtropical mountainous and hilly areas, broadleaved forests can grow longer under the same similar climatic conditions compared to coniferous forests [48]. The study also found a significantly stronger trend in the lengthening of growing season for trees than for shrubs and herbs, which is the same as the findings of Zhu et al. [49] but in contrast to those of Ge et al. [50], who reported that the interannual variation trend for trees in China from the 1960s to the 2000s was significantly weaker than for herbaceous plants, and this difference in trend was due to differences in the study area and the species of the plants themselves.

4.2. Asymmetry in Contributions of SOS and EOS Trends to LOS

We found the asymmetry in contributions of the SOS and EOS trends to LOS variations by counting the percentage of pixels with positive and negative C-index values. The results show that SOS trends control 48.4% of LOS variations and EOS trends control 51.4% of LOS variations, which shows a stronger association between EOS trends and LOS variations compared to SOS (Figure S1). Previous studies illustrated that the lengthened growing season was mainly driven by delayed autumn phenology, which is consistent with our results [14,38,49]. However, other researchers found that it is the changes in SOS, and not EOS, that dominate the changes in growing season length [51,52]. It can be seen that there are differences in previous studies regarding the attribution of LOS variations. To investigate the reasons for such differences, we divided the trends of SOS, EOS, and LOS into six change patterns (Table 2 and Figure 7). Our results show that in addition to the main growth pattern of Type I (SOS advanced, EOS delayed, and LOS lengthened), 15.2% of the regions had Type V (SOS delayed, EOS delayed, and LOS lengthened), and the delayed EOS was significant in this pattern. Another 12% of the regions showed Type II (SOS advanced, EOS advanced, and LOS lengthened) growth pattern, and advanced SOS was significant. However, since the percentage of the region of the growth pattern Type II is smaller than that of Types I and V, it is still the trend of delayed EOS that dominates the variation in LOS for the whole study area, leading to asymmetry of the relative contribution of SOS and EOS to LOS.
In addition, we found that ENF, EBF, MF, SL, and CL were all controlled by EOS trends, while the variations in LOS for two vegetation types, DBF and GL, were controlled by SOS trends (Table 5). As Figure 4e shows, the growing season lengths of DBF and GL are short, and there are previous studies demonstrating that the effect of EOS shifts on vegetation with short growing season cycles is insignificant [53]. The percentage of growth pattern type II is higher than other vegetation types in DBF and GL, and the advance in SOS is also significant, resulting in SOS dominating the variation in LOS (Table 4). Therefore, we suggest that the asymmetry in SOS and EOS trends contributing to LOS is related to vegetation types, and that future studies should focus on vegetation types to accurately model and predict vegetation phenology periods.

4.3. Analysis of the Drivers of Interannual Variations in LSP

Previous studies showed that the interaction of meteorological, soil, and biological factors influenced the interannual variability of LSP [6,54]. Our results suggest that SWP is the most important driver of interannual variations in SOS and EOS across the QMs (Figure 8). This is mainly because shortwave radiation compensates for the lack of chilling demand during plant physiological dormancy through day length, i.e., longer daylight hours, and has an critical effect on SOS by delaying the accumulation of abscisic acid and slowing down the rate of leaf senescence, and also on EOS [21]. We also found that SWP and MD contributed a total of 42.3% to the interannual variations in SOS, and besides SWP, MD was also an important driver of the interannual SOS variation (Table 6 and Table S1). This is due to the fluctuation of the time interval between SOS and MD, which depends on the different developmental stages of the phenology events to a large extent and on the specific differences in the life history of the plant, and needs to be explained by the phenotypic plasticity of the individual and the adaptation to the environment [55]. Therefore, we suggest that the effect of MD on interannual SOS variation varies considerably among vegetation individuals. STP also contributed 12.0% importance in explaining the interannual SOS variation (Table S1). This is mainly due to the increased soil temperature, which accelerated the rate of leaf tip emergence and whole leaf expansion, thus promoting SOS [56]. Moreover, SWP and PP together explain the importance of 43.2% of the interannual EOS variation (Table S1). This is mainly because preseason shortwave radiation and precipitation control the availability of sunlight and water in vegetation, respectively, and reduced precipitation affects water transport capacity, which limits the photosynthetic rate of leaf, leading to lower utilization of light and water conditions by plants and affecting the interannual variation in EOS [20,57]. The combined contribution of MD, SWP, and PP to interannual EOS variation was also found to be as high as 54.0% (Table S1), suggesting that the lifecycle of vegetation is strongly regulated by its own rhythms under improved hydrothermal conditions, and that biological rhythms play a critical role in interannual EOS variation [7].
Furthermore, our study shows that the effect of each driver on interannual variations in LSP was varied for different vegetation types (Figure 9 and Table 6). For example, SWP is the most important driver for ENF, MF, GL, and CL; MD is the most important driver for EBF and DBF; and STP is the most important driver for SL (Table 6). This difference is mainly due to the diversity of plant physiological structures and the different adaptive strategies of plants to environmental changes [58]. SWP has the greatest effect on the interannual SOS variation in GL, which is mainly distributed in higher parts of the QMs and receives abundant solar radiation. The strong solar radiation promotes root activity and advances SOS [59]. The greatest contribution of MD to interannual SOS variation in EBF is related to the fact that EBF grows mainly in the southern part of the QMs, where its deeply rooted system and water conservation adaptations combine to reduce water stress under the influence of a humid northern subtropical climate. This adaptation to environmental changes is strongly regulated by its own rhythms, such that EBF is most affected by MD [60]. The effect of STP in SL is mainly due to the preseason accumulation of soil temperatures susceptible to specific thresholds that accelerate soil thaw and vegetation wake, triggering SOS [61]. SWP was the main driver of interannual EOS variation for all vegetation types, with GL being most influenced by SWP (Table 6). This is because abundant solar radiation increases surface evaporation and reduces water availability in grasslands, which subsequently inhibits vegetation growth, resulting in the EOS of GL being most influenced by SWP [62]. PP and MD have the strongest effect on interannual EOS variation in SL (Table 6). To our best knowledge, SL is mainly distributed in semi-humid and semiarid areas, and the control of plant metabolism by water stress affects its transpiration and photosynthesis, resulting in impaired ATPase synthesis and accelerated chlorophyll degradation. Meanwhile the adaptation of vegetation to such adversity changes also affects the interannual variations in EOS [63]. The relationship between regional climate and vegetation phenology growth will be further explored in future studies.

4.4. Evaluation of RF Model

We validated the accuracy of our RF model (Figure S2). The R2 values of the linear regression formed by the predicted values and the observations inversions are both greater than 0.900 and 0.911, respectively. When predicting the dates of SOS or EOS, both RMSE and MAE are relatively small, showing that our RF model has good predictive performance. Machine learning, as a nonparametric multivariate approach, can integrate complex relationships between multiple spatial and temporal LSP dates and climate into a single model for predicting SOS and EOS. Then, the main drivers of SOS and EOS can be identified by estimating the importance of each variable [6]. However, it should be noted that although we have tried our best to adjust the hyperparameters of the algorithm to prevent overfitting in the model, some errors still appear in the test datasets (Table S2). Therefore, to obtain a better fit, it is necessary to further refine the study area and tree species in the future. Further study should compare different algorithms to better simulate the phenology period.

5. Conclusions

This study used the phenology metrics of vegetation in the QMs extracted from satellite NDVI data to analyze the spatiotemporal trends of LSP during 2001–2019, and to identify the dominant growth patterns of different vegetation types during the growing season. Furthermore, driving factors influencing interannual variations in LSP were emphatically investigated using the RF model. The main conclusions were as follows:
(1)
The average advance of SOS across QMs was 1.5 days/decade, with a significant advance in SOS observed for 27.5% of pixels. EOS was delayed by 2.4 days/decade, with a significant delay in EOS observed for 42.1% of pixels. LOS was lengthened by 3.9 days/decade, with a significant LOS lengthening observed for 40.3% of pixels.
(2)
The dominant pattern of change in the growing season for different vegetation types was advanced SOS, delayed EOS, and lengthened LOS, and this pattern had the highest percentage in evergreen broadleaved forests. The percentage of area shows that the patterns of delayed SOS and EOS and lengthened LOS were the highest percentage in shrubs.
(3)
For the whole QMs, LOS changes were mainly controlled by EOS, and the percentage was 51.4%. For deciduous broadleaved forests and grasses, LOS changes were attributed to SOS, while for other vegetation types, they were attributed to EOS.
(4)
SWP was found to be the most important factor influencing SOS and EOS, and grass and crop most influenced by SWP. Interannual variations in SOS were more influenced by biological factors (MD) than in EOS. The interannual variability of EOS is more influenced by preseason precipitation (PP) than SOS.

Supplementary Materials

The following are available online at https://www.mdpi.com/article/10.3390/rs13224538/s1, Figure S1: The primary factor (SOS or EOS) controlling change in LOS across the entire study area calculated by the C-index, Figure S2: Red fitting curves between predicted LSP dates from RF model and observed LSP dates from MODIS13A2-NDVI, Table S1: The relative importance of each driver affecting interannual variation in LSP, Table S2: Correlation results of random forest models built for different vegetation type cover areas.

Author Contributions

J.G.: Formal analysis, Methodology, Writing—original draft. X.L.: Conceptualization, Funding acquisition, Methodology, Supervision, Writing—original draft. W.G.: Data curation, Visualization. X.N.: Data curation, Visualization. W.M.: Writing—review & editing. Q.L.: Writing—review & editing. X.X.: Resources. All authors have read and agreed to the published version of the manuscript.

Funding

This research was supported by the National Natural Science Foundation of China (Grants 41971104) and partially supported by the Fundamental Research Funds for the Central Universities (Projects GK202107009 and 2021TS015).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The datasets generated and/or analyzed during the study are available from the corresponding author upon reasonable request.

Acknowledgments

I am very grateful to the reviewers and editors for their suggestions on the revision of the manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

LSPLand surface phenology
QMsQinling Mountains
NDVINormalized difference vegetation index
SOSThe start of the growing season
EOSThe end of the growing season
LOSThe length of the growing season
ENFEvergreen needleleaved forest
EBFEvergreen broadleaved forest
DBFDeciduous broadleaved forest
MFMixed forest
SLShrubland
GLGrassland
CLCropland
TPPreseason average temperature
TGGrowing season average temperature
PPPreseason total precipitation
PGGrowing season total precipitation
SWPPreseason mean shortwave radiation
SWGGrowing season mean shortwave radiation
STPPreseason soil temperature
STGGrowing season soil temperature
SMPPreseason soil moisture
SMGGrowing season soil moisture
MDThe middle date of the growing season
MNMaximum NDVI during growing season
RFRandom forest
OOBOut of bag

References

  1. Güsewell, S.; Furrer, R.; Gehrig, R.; Pietragalla, B. Changes in temperature sensitivity of spring phenology with recent climate warming in Switzerland are related to shifts of the preseason. Glob. Chang. Biol. 2017, 23, 5189–5202. [Google Scholar] [CrossRef]
  2. Fu, Y.H.; Piao, S.; Vitasse, Y.; Zhao, H.; De Boeck, H.J.; Liu, Q.; Yang, H.; Weber, U.; Hänninen, H.; Janssens, I.A. Increased heat requirement for leaf flushing in temperate woody species over 1980–2012: Effects of chilling, precipitation and insolation. Glob. Chang. Biol. 2015, 21, 2687–2697. [Google Scholar] [CrossRef]
  3. Lang, W.; Chen, X.; Qian, S.; Liu, G.; Piao, S. A new process-based model for predicting autumn phenology: How is leaf senescence controlled by photoperiod and temperature coupling? Agric. For. Meteorol. 2019, 268, 124–135. [Google Scholar] [CrossRef]
  4. Thackeray, S.J.; Henrys, P.A.; Hemming, D.; Bell, J.R.; Botham, M.S.; Burthe, S.; Helaouet, P.; Johns, D.G.; Jones, I.D.; Leech, D.I.; et al. Phenological sensitivity to climate across taxa and trophic levels. Nature 2016, 535, 241–245. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  5. Zeng, Z.; Piao, S.; Li, L.; Zhou, L.; Ciais, P.; Wang, T.; Li, Y.; Lian, X.; Wood, E.; Mao, J. Climate mitigation from vegetation biophysical feedbacks during the past three decades. Nat. Clim. Chang. 2017, 7, 432–436. [Google Scholar] [CrossRef]
  6. Rodriguez-Galiano, V.F.; Sanchez-Castillo, M.; Dash, J.; Atkinson, P.M.; Ojeda-Zujar, J. Modelling interannual variation in the spring and autumn land surface phenology of the European forest. Biogeosciences 2016, 13, 3305–3317. [Google Scholar] [CrossRef] [Green Version]
  7. Zu, J.; Zhang, Y.; Huang, K.; Liu, Y.; Chen, N.; Cong, N. Biological and climate factors co-regulated spatial-temporal dynamics of vegetation autumn phenology on the Tibetan Plateau. Int. J. Appl. Earth Obs. Geoinf. 2018, 69, 198–205. [Google Scholar] [CrossRef]
  8. Keenan, T.F.; Gray, J.; Friedl, M.A.; Toomey, M.; Bohrer, G.; Hollinger, D.Y.; Munger, J.W.; Keefe, J.O.; Schmid, H.P.; Wing, I.S. Net carbon uptake has increased through warming-induced changes in temperate forest phenology. Nat. Clim. Chang. 2014, 4, 598–604. [Google Scholar] [CrossRef]
  9. Wu, C.; Hou, X.; Peng, D.; Gonsamo, A.; Xu, S. Land surface phenology of China’s temperate ecosystems over 1999–2013: Spatial–temporal patterns, interaction effects, covariation with climate and implications for productivity. Agric. For. Meteorol. 2016, 216, 177–187. [Google Scholar] [CrossRef]
  10. Beurs, K.M.D.; Henebry, G.M. Land surface phenology and temperature variation in the International Geosphere-Biosphere Program high-latitude transects. Glob. Chang. Biol. 2010, 11, 779–790. [Google Scholar] [CrossRef]
  11. Cleland, E.E.; Chuine, I.; Menzel, A.; Mooney, H.A.; Schwartz, M.D. Shifting plant phenology in response to global change. Trends Ecol. Evol. 2007, 22, 357–365. [Google Scholar] [CrossRef]
  12. Stoeckli, R.; Vidale, P.L. European plant phenology and climate as seen in a 20-year AVHRR land-surface parameter dataset. Int. J. Remote Sens. 2004, 25, 3303–3330. [Google Scholar] [CrossRef]
  13. Piao, S.; Fang, J.; Zhou, L.; Ciais, P.; Zhu, B. Variations in satellite-derived phenology in China’s temperate vegetation. Glob. Chang. Biol. 2010, 12, 672–685. [Google Scholar] [CrossRef]
  14. Jeong, S.J.; Chang-Hoi, H.O.; Gim, H.J.; Brown, M.E. Phenology shifts at start vs. end of growing season in temperate vegetation over the Northern Hemisphere for the period 1982–2008. Glob. Chang. Biol. 2011, 17, 2385–2399. [Google Scholar] [CrossRef]
  15. Shen, M.; Tang, Y.; Chen, J.; Zhu, X.; Zheng, Y. Influences of temperature and precipitation before the growing season on spring phenology in grasslands of the central and eastern Qinghai-Tibetan Plateau. Agric. For. Meteorol. 2011, 151, 1711–1722. [Google Scholar] [CrossRef]
  16. Yuan, M.; Zhao, L.; Lin, A.; Li, Q.; Qu, S. How do climatic and non-climatic factors contribute to the dynamics of vegetation autumn phenology in the Yellow River Basin, China? Ecol. Indic. 2020, 112, 106112. [Google Scholar] [CrossRef]
  17. Wang, C.; Guo, H.; Zhang, L.; Liu, S.; Qiu, Y.; Sun, Z. Assessing phenological change and climatic control of alpine grasslands in the Tibetan Plateau with MODIS time series. Int. J. Biometeorol. 2015, 59, 11–23. [Google Scholar] [CrossRef]
  18. Flynn, D.; Wolkovich, E.M. Temperature and photoperiod drive spring phenology across all species in a temperate forest community. New Phytol. 2018, 219, 1353–1362. [Google Scholar] [CrossRef] [Green Version]
  19. Liu, Q.; Piao, S.; Fu, Y.H.; Gao, M.; Peñuelas, J.; Janssens, I.A. Climatic Warming Increases Spatial Synchrony in Spring Vegetation Phenology Across the Northern Hemisphere. Geophys. Res. Lett. 2019, 46, 1641–1650. [Google Scholar] [CrossRef] [Green Version]
  20. Prevéy, J.S.; Seastedt, T.R. Seasonality of precipitation interacts with exotic species to alter composition and phenology of a semi-arid grassland. J. Ecol. 2014, 102, 1549–1561. [Google Scholar] [CrossRef]
  21. Way, D.A.; Montgomery, R.A. Photoperiod constraints on tree phenology, performance and migration in a warming world. Plant Cell Environ. 2015, 38, 1725–1736. [Google Scholar] [CrossRef] [PubMed]
  22. Ztürk, M.; Bolat, L.; Ergün, A. Influence of air-soil temperature on leaf expansion and LAI of Carpinus betulus trees in a temperate urban forest patch. Agric. For. Meteorol. 2015, 200, 185–191. [Google Scholar] [CrossRef]
  23. Jin, Z.; Zhuang, Q.; He, J.S.; Luo, T.; Shi, Y. Phenology shift from 1989 to 2008 on the Tibetan Plateau: An analysis with a process-based soil physical model and remote sensing data. Clim. Chang. 2013, 119, 435–449. [Google Scholar] [CrossRef]
  24. Zhou, Y. Asymmetric Behavior of Vegetation Seasonal Growth and the Climatic Cause: Evidence from Long-Term NDVI Dataset in Northeast China. Remote Sens. 2019, 11, 2107. [Google Scholar] [CrossRef] [Green Version]
  25. Yang, M.; Nelson, F.E.; Shiklomanov, N.I.; Guo, D.; Wan, G. Permafrost degradation and its environmental effects on the Tibetan Plateau: A review of recent research. Earth-Sci. Rev. 2010, 103, 31–44. [Google Scholar] [CrossRef]
  26. Park, H.; Jeong, S.J.; Ho, C.H.; Kim, J.; Brown, M.E.; Schaepman, M.E. Nonlinear response of vegetation green-up to local temperature variations in temperate and boreal forests in the Northern Hemisphere. Remote Sens. Environ. 2015, 165, 100–108. [Google Scholar] [CrossRef]
  27. Breiman, L. Random forests. Mach. Learn. 2001, 45, 5–32. [Google Scholar] [CrossRef] [Green Version]
  28. Li, Z.; Liang, M.; Li, Z.; Pierre, M.; Tong, X.; Zhang, J.; Dong, L.; Zheng, Y.; Ma, W.; Zhao, L. Plant functional groups mediate effects of climate and soil factors on species richness and community biomass in grasslands of Mongolian Plateau. J. Plant Ecol. 2021, 14, 679–691. [Google Scholar] [CrossRef]
  29. Zhao, M.; Peng, C.; Xiang, W.; Deng, X.; Tian, D.; Zhou, X.; Yu, G.; He, H.; Zhao, Z. Plant phenological modeling and its application in global climate change research: Overview and future challenges. Environ. Rev. 2013, 21, 1–14. [Google Scholar] [CrossRef]
  30. Deng, C.; Bai, H.; Ma, X.; Zhao, T.; Huang, X. Spatiotemporal differences in the climatic growing season in the Qinling Mountains of China under the influence of global warming from 1964 to 2015. Theor. Appl. Climatol. 2019, 138, 1899–1911. [Google Scholar] [CrossRef]
  31. Bai, W.G.; Zhang, L. Phytocoenological characteristics and community classification of Abies fargesii forests in Qinling Mountains. J. Beijing For. Univ. 2007, S2, 222–226. [Google Scholar]
  32. Myneni, R.B.; Keeling, C.D.; Tucker, C.J.; Asrar, G.; Nemani, R.R. Increased plant growth in the northern high latitudes from 1981 to 1991. Nature 1997, 386, 698–702. [Google Scholar] [CrossRef]
  33. Wang, X.; Xiao, J.; Li, X.; Cheng, G.; Ma, M.; Zhu, G.; Altaf Arain, M.; Andrew Black, T.; Jassal, R.S. No trends in spring and autumn phenology during the global warming hiatus. Nat. Commun. 2019, 10, 2389. [Google Scholar] [CrossRef] [PubMed]
  34. Cao, R.; Yang, C.; Shen, M.; Chen, J.; Ji, Z.; Cong, W.; Wei, Y. A simple method to improve the quality of NDVI time-series data by integrating spatiotemporal information with the Savitzky-Golay filter. Remote Sens. Environ. 2018, 217, 244–257. [Google Scholar] [CrossRef]
  35. Kaur, R.; Kiran, G.S.; Shah, M.N.; Mistry, N.V.; Mohan, S. Applicability of Smoothing Techniques in Generation of Phenological Metrics of Tectona grandis L. Using NDVI Time Series Data. Remote Sens. 2021, 13, 3343. [Google Scholar]
  36. Gocic, M.; Trajkovic, S. Analysis of changes in meteorological variables using Mann-Kendall and Sen’s slope estimator statistical tests in Serbia. Glob. Planet. Chang. 2013, 100, 172–182. [Google Scholar] [CrossRef]
  37. Mann, H.B. Nonparametric test against trend. Econometrica 1945, 13, 245–259. [Google Scholar] [CrossRef]
  38. Garonna, I.; Jong, R.D.; Wit, A.D.; Mücher, C.A.; Schmid, B.; Schaepman, M.E. Strong contribution of autumn phenology to changes in satellite-derived growing season length estimates across Europe (1982–2011). Glob. Chang. Biol. 2015, 20, 3457–3470. [Google Scholar] [CrossRef]
  39. Abdel-Rahman, E.M.; Ahmed, F.B.; Ismail, R. Random forest regression and spectral band selection for estimating sugarcane leaf nitrogen concentration using EO-1 Hyperion hyperspectral data. Int. J. Remote Sens. 2013, 34, 712–728. [Google Scholar] [CrossRef]
  40. Tian, F.; Cai, Z.; Jin, H.; Hufkens, K.; Scheifinger, H.; Tagesson, T.; Smets, B.; Van Hoolst, R.; Bonte, K.; Ivits, E.; et al. Calibrating vegetation phenology from Sentinel-2 using eddy covariance, PhenoCam, and PEP725 networks across Europe. Remote Sens. Environ. 2021, 260, 112456. [Google Scholar] [CrossRef]
  41. Sun, Q.; Li, B.; Zhou, G.; Jiang, Y.; Yuan, Y. Delayed autumn leaf senescence date prolongs the growing season length of herbaceous plants on the Qinghai-Tibetan Plateau. Agric. For. Meteorol. 2020, 284, 107896. [Google Scholar] [CrossRef]
  42. Wang, X.; Piao, S.; Xu, X.; Ciais, P.; Macbean, N.; Myneni, R.B.; Li, L. Has the advancing onset of spring vegetation green-up slowed down or changed abruptly over the last three decades? Glob. Ecol. Biogeogr. 2015, 24, 621–631. [Google Scholar] [CrossRef]
  43. Zeng, H.; Jia, G.; Epstein, H. Recent changes in phenology over the northern high latitudes detected from multi-satellite data. Environ. Res. Lett. 2011, 6, 45508–45518. [Google Scholar] [CrossRef] [Green Version]
  44. Xia, H.; Qin, Y.; Feng, G.; Meng, Q.; Liu, G. Forest Phenology Dynamics to Climate Change and Topography in a Geographic and Climate Transition Zone: The Qinling Mountains in Central China. Forests 2019, 10, 1007. [Google Scholar] [CrossRef] [Green Version]
  45. Wang, H.; Ge, Q.; Dai, J.; Tao, Z. Geographical pattern in first bloom variability and its relation to temperature sensitivity in the USA and China. Int. J. Biometeorol. 2015, 59, 961–969. [Google Scholar] [CrossRef]
  46. Panchen, Z.A.; Primack, R.B.; Nordt, B.; Ellwood, E.R.; Stevens, A.D.; Renner, S.S.; Willis, C.G.; Fahey, R.; Whittemore, A.; Du, Y. Leaf out times of temperate woody plants are related to phylogeny, deciduousness, growth habit and wood anatomy. New Phytol. 2014, 203, 1208–1219. [Google Scholar] [CrossRef]
  47. Kurukulasuriya, P.; Mendelsohn, R. Impact And Adaptation Of South-East Asian Farmers To Climate Change: Conclusions And Policy Recommendations. Clim. Chang. Econ. 2017, 8, 1740007. [Google Scholar] [CrossRef]
  48. Qiu, B.W.; Zhong, M. Spatiotemporal variability of vegetation phenology with reference to altitude and climate in the subtropical mountain and hill region, China. Chin. Sci. Bull. 2013, 58, 2883–2892. [Google Scholar] [CrossRef] [Green Version]
  49. Zheng, Z.; Zhu, W.; Chen, G.; Jiang, N.; Fan, D.; Zhang, D. Continuous but diverse advancement of spring-summer phenology in response to climate warming across the Qinghai-Tibetan Plateau. Agric. For. Meteorol. 2016, 223, 194–202. [Google Scholar] [CrossRef]
  50. Ge, Q.; Wang, H.; Rutishauser, T.; Dai, J. Phenological response to climate change in China: A meta-analysis. Glob. Chang. Biol. 2015, 21, 265–274. [Google Scholar] [CrossRef]
  51. Park, T.; Ganguly, S.; T Mmervik, H.; Euskirchen, E.S.; H Gda, K.A.; Karlsen, S.R.; Brovkin, V.; Nemani, R.R.; Myneni, R.B. Changes in growing season duration and productivity of northern vegetation inferred from long-term remote sensing data. Environ. Res. Lett. 2016, 11, 84001. [Google Scholar] [CrossRef]
  52. Shen, M.; Piao, S.; Dorji, T.; Liu, Q.; Cong, N.; Chen, X.; An, S.; Wang, S.; Wang, T.; Zhang, G. Plant phenological responses to climate change on the Tibetan Plateau: Research status and challenges. Natl. Sci. Rev. 2015, 2, 454–467. [Google Scholar] [CrossRef] [Green Version]
  53. Cong, N.; Huang, K.; Zhang, Y. Unsynchronized Driving Mechanisms of Spring and Autumn Phenology Over Northern Hemisphere Grasslands. Front. For. Glob. Chang. 2021, 3, 610162. [Google Scholar] [CrossRef]
  54. Wang, J.; Sun, H.; Xiong, J.; He, D.; Cheng, W.; Ye, C.; Yong, Z.; Huang, X. Dynamics and Drivers of Vegetation Phenology in Three-River Headwaters Region Based on the Google Earth Engine. Remote Sens. 2021, 13, 2528. [Google Scholar] [CrossRef]
  55. Mcwatters, H.G.; Devlin, P.F. Timing in plants—A rhythmic arrangement. Febs. Lett. 2011, 585, 1474–1484. [Google Scholar] [CrossRef] [Green Version]
  56. Stone, P.J.; Sorensen, I.B.; Jamieson, P.D. Effect of soil temperature on phenology, canopy development, biomass and yield of maize in a cool-temperate climate. Field Crop. Res. 1999, 63, 169–178. [Google Scholar] [CrossRef]
  57. Chen, X.; Ciais, P.; Maignan, F.; Zhang, Y.; Zhang, H. Vapor Pressure Deficit and Sunlight Explain Seasonality of Leaf Phenology and Photosynthesis Across Amazonian Evergreen Broadleaved Forest. Glob. Biogeochem. Cycles 2021, 35, e2020GB006893. [Google Scholar] [CrossRef]
  58. Diez, J.M.; Ibáñez, I.; Miller-Rushing, A.J.; Mazer, S.J.; Crimmins, T.M.; Crimmins, M.A.; Bertelsen, C.D.; Inouye, D.W. Forecasting phenology: From species variability to community patterns. Ecol. Lett. 2012, 15, 545–553. [Google Scholar] [CrossRef] [Green Version]
  59. Li, K.; Wang, C.; Sun, Q.; Rong, G.; Tong, Z.; Liu, X.; Zhang, J. Spring Phenological Sensitivity to Climate Change in the Northern Hemisphere: Comprehensive Evaluation and Driving Force Analysis. Remote Sens. 2021, 10, 1972. [Google Scholar] [CrossRef]
  60. Zhao, A.; Zhang, A.; Cao, S.; Liu, X.; Liu, J.; Cheng, D. Responses of vegetation productivity to multi-scale drought in Loess Plateau, China. Catena 2017, 163, 165–171. [Google Scholar] [CrossRef]
  61. Chuine, I. A unified model for budburst of trees. J. Theor. Biol. 2000, 207, 337–347. [Google Scholar] [CrossRef] [PubMed]
  62. Chen, B.Z.; Zhang, H.F.; Yan, J.W.; Wang, G.Y.; Zhou, T.M.; Xu, G.; Innes, J.L.; Che, M.L.; Dou, X.M. Spatial and temporal variations in the end date of the vegetation growing season throughout the Qinghai-Tibetan Plateau from 1982 to 2011. Agric. For. Meteorol. 2014, 189–190, 81–90. [Google Scholar]
  63. Gerten, D.; Schaphoff, S.; Lucht, W. Potential future changes in water limitations of the terrestrial biosphere. Clim. Chang. 2007, 80, 277–299. [Google Scholar] [CrossRef]
Remotesensing 13 04538 g001 550
Figure 1. An overview of (a) the study area and (b) spatial distribution of vegetation types. Only pixels with unchanged vegetation types were analyzed in this study period. The dataset with a total of 218,899 pixels was divided into seven vegetation types.
Figure 1. An overview of (a) the study area and (b) spatial distribution of vegetation types. Only pixels with unchanged vegetation types were analyzed in this study period. The dataset with a total of 218,899 pixels was divided into seven vegetation types.
Remotesensing 13 04538 g001
Remotesensing 13 04538 g002 550
Figure 2. The description of phenology metrics correlations extracted using the NDVI time series datasets.
Figure 2. The description of phenology metrics correlations extracted using the NDVI time series datasets.
Remotesensing 13 04538 g002
Remotesensing 13 04538 g003 550
Figure 3. Flowchart of the method used in this study. MODIS 13A2 NDVI (2001–2019), in which the spatiotemporal distributions of NDVI for 2001, 2010, 2015 and 2019 are presented as examples of NDVI time series used as extract phenology metrics. Land cover (2001–2019), in which the spatiotemporal distribution of land cover for 2001, 2010, 2015 and 2019, used as an example of the changes in different vegetation types. Environmental factors, the first column indicates the drivers for SOS and the second column indicates the drivers for EOS, and shows the spatial distribution of each driver separately, which is used to example the predictors used for RF models.
Figure 3. Flowchart of the method used in this study. MODIS 13A2 NDVI (2001–2019), in which the spatiotemporal distributions of NDVI for 2001, 2010, 2015 and 2019 are presented as examples of NDVI time series used as extract phenology metrics. Land cover (2001–2019), in which the spatiotemporal distribution of land cover for 2001, 2010, 2015 and 2019, used as an example of the changes in different vegetation types. Environmental factors, the first column indicates the drivers for SOS and the second column indicates the drivers for EOS, and shows the spatial distribution of each driver separately, which is used to example the predictors used for RF models.
Remotesensing 13 04538 g003
Remotesensing 13 04538 g004 550
Figure 4. (a,c,e) Spatial distribution of the average phenology metrics from 2001 to 2019 and (b,d,f) standard deviation (Std) of the phenology metrics. Insets at bottom left show the histogram of the average pixel values for different vegetation types.
Figure 4. (a,c,e) Spatial distribution of the average phenology metrics from 2001 to 2019 and (b,d,f) standard deviation (Std) of the phenology metrics. Insets at bottom left show the histogram of the average pixel values for different vegetation types.
Remotesensing 13 04538 g004
Remotesensing 13 04538 g005 550
Figure 5. (a,c,e) Spatial distribution of phenology metrics trends from 2001 to 2019 and (b,d,f) significant (p < 0.1) changes in phenology metrics trends for the study periods. Insets at bottom left show the average pixel values of the trends of the phenology metrics for the different vegetation types.
Figure 5. (a,c,e) Spatial distribution of phenology metrics trends from 2001 to 2019 and (b,d,f) significant (p < 0.1) changes in phenology metrics trends for the study periods. Insets at bottom left show the average pixel values of the trends of the phenology metrics for the different vegetation types.
Remotesensing 13 04538 g005
Remotesensing 13 04538 g006 550
Figure 6. Interannual variations and significance of the mean phenology metrics for the entire study area and the areas covered by different vegetation types. The unit of Sen’s slope is days/year.
Figure 6. Interannual variations and significance of the mean phenology metrics for the entire study area and the areas covered by different vegetation types. The unit of Sen’s slope is days/year.
Remotesensing 13 04538 g006
Remotesensing 13 04538 g007 550
Figure 7. Trends and significance of each pattern during the growing season for different vegetation types. (a) growing season pattern I. (b) growing season pattern II. (c) growing season pattern III. (d) growing season pattern IV. (e) growing season pattern V. (f) growing season pattern VI. Note: “All” refers to all vegetation in the QMs.
Figure 7. Trends and significance of each pattern during the growing season for different vegetation types. (a) growing season pattern I. (b) growing season pattern II. (c) growing season pattern III. (d) growing season pattern IV. (e) growing season pattern V. (f) growing season pattern VI. Note: “All” refers to all vegetation in the QMs.
Remotesensing 13 04538 g007
Remotesensing 13 04538 g008 550
Figure 8. Relative importance of drivers of affecting interannual variations in LSP for the entire study area, in decreasing order. (a) The ranking of the drivers of affecting interannual variations in SOS. (b) The ranking of the drivers of affecting interannual variations in EOS. The specific values for the relative importance of each driver are in Table S1. Note: The abbreviated variable names are the same as in Table 3.
Figure 8. Relative importance of drivers of affecting interannual variations in LSP for the entire study area, in decreasing order. (a) The ranking of the drivers of affecting interannual variations in SOS. (b) The ranking of the drivers of affecting interannual variations in EOS. The specific values for the relative importance of each driver are in Table S1. Note: The abbreviated variable names are the same as in Table 3.
Remotesensing 13 04538 g008
Remotesensing 13 04538 g009 550
Figure 9. The relative importance of drivers that affect interannual variations in LSP for different vegetation types. Relative importance was derived from the importance scores of variables (VI) based on RF models established for different vegetation types: (a) SOS. (b) EOS. Different colors indicate different factors. Note: SOS or EOS as response variable in the RF model. The specific values for the relative importance of each driver are in Table S1.
Figure 9. The relative importance of drivers that affect interannual variations in LSP for different vegetation types. Relative importance was derived from the importance scores of variables (VI) based on RF models established for different vegetation types: (a) SOS. (b) EOS. Different colors indicate different factors. Note: SOS or EOS as response variable in the RF model. The specific values for the relative importance of each driver are in Table S1.
Remotesensing 13 04538 g009
Table
Table 1. Datasets and sources in the study area.
Table 1. Datasets and sources in the study area.
DatasetSpatial ResolutionTemporal ResolutionTime SpanSource
MODIS13A2 NDVI1 km16 days2001–2019The Level-1 and Atmosphere Archive and Distribution System Distributed Active Archive Center (LAADS DAAC) (https://search.earthdata.nasa.gov/search/, accessed on 15 October 2020).
Land cover
(CCI-LC)		300 mYearly2001–2019http://maps.elie.ucl.ac.be/CCI/viewer/index.php, accessed on 20 October 2020
Temperature0.1°hourly2001–2019The Reanalysis (ERA5) climatic datasets (https://cds.climate.copernicus.eu, accessed on 12 November 2020)
Precipitation0.1°hourly2001–2019
Shortwave radiation0.1°hourly2001–2019
Soil temperature0.1°hourly2001–2019
Soil moisture0.1°hourly2001–2019
Table
Table 2. Six change patterns of plant growing seasons, i.e., six combination types of SOS, EOS, and LOS change trends. The plus and minus signs represent the trend direction corresponding to SOS, EOS, or LOS, respectively.
Table 2. Six change patterns of plant growing seasons, i.e., six combination types of SOS, EOS, and LOS change trends. The plus and minus signs represent the trend direction corresponding to SOS, EOS, or LOS, respectively.
Change PatternTrend of SOSTrend of EOSTrend of LOS
IAdvanced (−)Delayed (+)Lengthened (+)
IIAdvanced (−)Advanced (−)Lengthened (+)
IIIAdvanced (−)Advanced (−)Shortened (−)
IVDelayed (+)Advanced (−)Shortened (−)
VDelayed (+)Delayed (+)Lengthened (+)
VIDelayed (+)Delayed (+)Shortened (−)
Table
Table 3. Predictive variables used in the modeling of the LSP dates. The 12 predictive variables were classified into three categories: meteorological factors, soil factors, and biological factors. Meteorological factors include TP, TG, PP, PG, SWP, and SWG (6 in total). Soil factors include STP, STG, SMP, and SMG (4 in total). Biological factors include MN and MD (2 in total). Growing season is defined as the days between SOS and EOS.
Table 3. Predictive variables used in the modeling of the LSP dates. The 12 predictive variables were classified into three categories: meteorological factors, soil factors, and biological factors. Meteorological factors include TP, TG, PP, PG, SWP, and SWG (6 in total). Soil factors include STP, STG, SMP, and SMG (4 in total). Biological factors include MN and MD (2 in total). Growing season is defined as the days between SOS and EOS.
VariablesSOS Drivers EOS Drivers
Meteorological factorsPreseason average temperature * (TP)Preseason average temperature ** (TP)
Growing season average temperature (TG)Growing season average temperature (TG)
Preseason total precipitation * (PP)Preseason total precipitation ** (PP)
Growing season total precipitation (PG)Growing season total precipitation (PG)
Preseason mean shortwave radiation * (SWP)Preseason mean shortwave radiation ** (SWP)
Growing season mean shortwave radiation (SWG)Growing season mean shortwave radiation (SWG)
Soil factorsPreseason soil temperature * (STP)Preseason soil temperature ** (STP)
Growing season soil temperature (STG)Growing season soil temperature (STG)
Preseason soil moisture * (SMP)Preseason soil moisture ** (SMP)
Growing season soil moisture (SMG)Growing season soil moisture (SMG)
Biological factorsMaximum NDVI during growing season (MN)Maximum NDVI during growing season (MN)
Middle season date (MD)Middle season date (MD)
* Predicted over a period of 30 days. ** Predicted over a period of 60 days.
Table
Table 4. The percentage of phenology metrics datasets consisting of each pixel showing different change patterns in the growing seasons for all the vegetation types. Change patterns refer to the trend groupings in Table 2.
Table 4. The percentage of phenology metrics datasets consisting of each pixel showing different change patterns in the growing seasons for all the vegetation types. Change patterns refer to the trend groupings in Table 2.
Change PatternsAll the Vegetation TypesENFEBFDBFMFSLGLCL
I48.4%46.6%50.1%49.3%47.2%43.2%44.8%48.7%
II12.0%10.9%9.5%13.7%11.8%6.9%15.2%10.5%
III7.3%6.8%6.2%8.6%7.3%3.1%7.6%5.9%
IV8.7%9.0%7.7%7.9%8.8%15.5%12.2%9.3%
V15.2%16.8%18.0%12.4%16.5%23.5%12.7%17.2%
VI8.5%9.9%8.5%8.1%8.4%7.8%7.5%8.6%
Total100.0%100.0%100.0%100.0%100.0%100.0%100.0%100.0%
Table
Table 5. The percentage of datasets in which LOS change was primarily attributable to the shift in SOS or EOS. For most vegetation types, the percentages show that the trend in LOS was mainly controlled by the shift in EOS.
Table 5. The percentage of datasets in which LOS change was primarily attributable to the shift in SOS or EOS. For most vegetation types, the percentages show that the trend in LOS was mainly controlled by the shift in EOS.
All the Vegetation TypesSOS ControlledEOS ControlledTotal
ENF46.9%53.1%100%
EBF43.4%56.6%100%
DBF53.4%46.6%100%
MF47.2%52.8%100%
SL41.7%58.3%100%
GL52.0%48.0%100%
CL45.9%54.1%100%
The whole area48.6%51.4%100%
Table
Table 6. The top three dominant drivers affecting interannual variations in LSP. These three dominant factors are derived from the ranking of the importance scores of the variables (VI). Different vegetation types have different dominant drivers.
Table 6. The top three dominant drivers affecting interannual variations in LSP. These three dominant factors are derived from the ranking of the importance scores of the variables (VI). Different vegetation types have different dominant drivers.
LSPAll the Vegetation TypesFirst Dominant DriverSecond Dominant DriverThird Dominant Driver
SOSENFSWPMDSTP
EBFMDSWPPP
DBFMDSWPSTP
MFSWPMDSTP
SLSTPTPSWP
GLSWPMDSTP
CLSWPMDTG
the whole areaSWPMDSTP
EOSENFSWPTPPP
EBFSWPMDPP
DBFSWPPPTP
MFSWPPPMD
SLSWPPPMD
GLSWPSWGTP
CLSWPTPPP
the whole areaSWPPPMD
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Guo, J.; Liu, X.; Ge, W.; Ni, X.; Ma, W.; Lu, Q.; Xing, X. Specific Drivers and Responses to Land Surface Phenology of Different Vegetation Types in the Qinling Mountains, Central China. Remote Sens. 2021, 13, 4538. https://doi.org/10.3390/rs13224538

AMA Style

Guo J, Liu X, Ge W, Ni X, Ma W, Lu Q, Xing X. Specific Drivers and Responses to Land Surface Phenology of Different Vegetation Types in the Qinling Mountains, Central China. Remote Sensing. 2021; 13(22):4538. https://doi.org/10.3390/rs13224538

Chicago/Turabian Style

Guo, Jiaqi, Xiaohong Liu, Wensen Ge, Xiaofeng Ni, Wenyuan Ma, Qiangqiang Lu, and Xiaoyu Xing. 2021. "Specific Drivers and Responses to Land Surface Phenology of Different Vegetation Types in the Qinling Mountains, Central China" Remote Sensing 13, no. 22: 4538. https://doi.org/10.3390/rs13224538

APA Style

Guo, J., Liu, X., Ge, W., Ni, X., Ma, W., Lu, Q., & Xing, X. (2021). Specific Drivers and Responses to Land Surface Phenology of Different Vegetation Types in the Qinling Mountains, Central China. Remote Sensing, 13(22), 4538. https://doi.org/10.3390/rs13224538

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop