**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).

**Figure 4.** (**<sup>a</sup>**,**c**,**<sup>e</sup>**) 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.

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).

**Figure 5.** (**<sup>a</sup>**,**c**,**<sup>e</sup>**) 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.

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).

**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.

#### *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.

**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.


**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.

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).

**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.


*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%).

**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.** *Cont.*

**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.

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.

**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.

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.
