*2.4. Methods*

2.4.1. Extraction of Plant Phenological Information

(1) *NDVI*

In this study, on the basis of *NDVI* that is estimated by the MOD09A1 band information, we calculated the SOS and EOS using the method of the relative and absolute rates of *NDVI* change, respectively. The *NDVI* is defined by [31]:

$$NDVI = (\rho\_{NIR} - \rho\_{Red}) / (\rho\_{NIR} + \rho\_{Red}),\tag{1}$$

where *ρNIR* and *ρRed* are the spectral reflectance values calculated in the near-infrared and red bands, respectively.

(2) Determination of the SOS and EOS

We used the maximum relative and minimum absolute rates of change in *NDVI* to calculate the SOS and EOS based on previous studies [14]. The equations of these rates of change can be expressed by:

$$NDVI\_{rate\\_rel} = \frac{NDVI\_{t+1} - NDVI\_{t}}{NDVI\_{t}}, \ t \in [1, 2, \dots, 365],\tag{2}$$

$$NDVI\_{\text{rate\\_abs}} = NDVI\_{t+1} - NDVI\_{t\\_t} \; t \in [1, 2, \dots, 365],\tag{3}$$

where *NDVIrate\_rel* and *NDVIrate\_abs* are the relative and absolute rates of change, respectively.

The specific calculation process is as follows. First, we calculated the time (*T*) when the maximum value appears based on the *NDVI* time-series data. The *NDVI* curve was divided into a rising (0, T) and a descending (T, 365) stage. Second, based on Equations (2) and (3), the maximum relative and minimum absolute rates of change were calculated by using the *NDVI* time-series data. Then, the thresholds of SOS and EOS were determined based on the maximum relative and minimum absolute rates of change, respectively. Third, if the *NDVI* value was greater than the SOS threshold at time 0 to *T*, the corresponding date of the year was regarded as the SOS. Similarly, if the *NDVI* value of some pixels was less than the EOS threshold at time *T* to 365, the corresponding day of the year plus one was regarded as the EOS (Figure 3b).

#### 2.4.2. The Spatiotemporal Pattern of Plant Phenology

#### (1) Linear Regression Analysis

We adopted a linear regression analysis to analyze the monotonic trend of the vegetation phenology and indicators [32,33]. The trend slope in a multi-year regression equation represents the amount of inter-annual change and can be found using the least squares method as follows:

$$Slope = \frac{n \cdot \sum\_{t=1}^{n} t \cdot X\_t - \sum\_{t=1}^{n} t \cdot \sum\_{t=1}^{n} X\_t}{n \cdot \sum\_{t=1}^{n} t^2 - \left(\sum\_{t=1}^{n} t\right)^2},\tag{4}$$

where *Slope* refers to the inter-annual trend, *n* is the number of years of the study, and the *Xt* is the value of this variable in the *t*-th year. When the slope is positive or negative, this indicates an increasing or decreasing trend, respectively.

#### (2) Standard Deviation Analysis

Standard deviation is a measure of the degree of data dispersion that can reflect the stability or fluctuation of variables [34]. For this study, the stability or fluctuation of plant phenology was calculated by standard deviation based on the pixel scale. The calculation formula is as follows:

$$S\_i = \sqrt{\frac{1}{n} \sum\_{i=1}^{n} \left(X\_i - \overline{X}\right)^2},\tag{5}$$

where *Si* indicates the standard deviation of an *X* dataset. When the *Si* value is larger, the distribution of the data is more discrete and has a larger range of fluctuation. In contrast, when the *Si* value is smaller, the distribution of the data is more concentrated and the range of fluctuation is smaller.

#### 2.4.3. Driving Force Analysis

(1) Pearson Correlation Coefficient

For this paper, we used correlation analysis to determine the relationship between the plant phenology (SOS, EOS, and LOS) and other factors. A higher value indicates a stronger correlation; otherwise, it means a weaker correlation [28,35]. The relevant formula is as follows:

$$R\_{xy} = \frac{\sum\_{i=1}^{n} \left[ \left[ \mathbf{x}\_i - \overline{\mathbf{x}} \right] \cdot \left[ \mathbf{y}\_i - \overline{\mathbf{y}} \right] \right]}{\sqrt{\sum\_{i=1}^{n} \left[ \left[ \mathbf{x}\_i - \overline{\mathbf{x}} \right]^2 \cdot \left[ \mathbf{y}\_i - \overline{\mathbf{y}} \right]^2 \right]}},\tag{6}$$

where *Rxy* is the correlation coefficient between *x* and *y*, *n* is the number of years during the study, *xi* and *yi* are the two sets of variables, and *x* and *y* are the mean values of variables.

#### (2) Structural Equation Model

SEM is a method used to analyze the relationship between variables based on a covariance matrix of variables, which includes maximum likelihood, synthesis of factor, and path analyses [24]. It pre-sets the dependence relationship between the factors in the system based on the researcher's prior knowledge, which can judge the strength of the relationship between the factors and can fit and judge the overall model. In addition, SEM has several advantages. For example, the direct or indirect effects of a particular variable on another variable can be partitioned by SEM, and SEM estimates and reports the total

path coefficient to present the strengths of these multiple effects [36]. Since the change in SOS and EOS eventually lead to the change in LOS, this paper only used SEM to explore the potential influence mechanism of climate and soil factors on LOS in the TRHR.
