Vegetation Phenology and Its Response to Climate Change in Extremely Arid Area: A Case Study of Populus Euphratica in The Upper Tarim River Basin, NW China

Hualin Li College of Soil and Water Conservation, Beijing Forestry University Jianzhong Feng Agricultural Information Institute, Chinese Academy of Agricultural Sciences Linyan bai (  baily@aircas.ac.cn ) Key Laboratory of Digital Earth Sciences, Aerospace Information Research Institute, Chinese Academy of Sciences https://orcid.org/0000-0002-9376-3255 Jianjun Zhang College of Soil and Water Conservation, Beijing Forestry University


Introduction
Vegetation phenology is the research of events in the plant life cycle, and how these events adapt to climate changes (Zhang et  The plant growth such as leaf unfolding, owering, fruiting, and leaf senescence is recorded to re ect the phenological characteristics of vegetation based on ground monitoring sites of phenology. However, this method has some limitations in spatiotemporal scale and biome-scale (Zhang et al., 2021a). Therefore, the remote sensing technique for phenological monitoring, which is used to simulate various logistics regression models through long-term series remote sensing data and to extract vegetation phenology (i.e., the start of the growing season (SOS), length of the growing season (LOS), and end of the growing season (EOS)) from suitable models, has been widely applied to explore the spatial-temporal variations of vegetation phenology at the regional and global scales (Steven P Arid and semi-arid ecosystems are important parts of the Earth's climate system, covering 41% of the Earth's land surface (Schimel, 2010;Yuan et al., 2014;Yu et al., 2021). These regions are characterized by harsh environments and little rainfall, especially in highly arid areas, where the average annual rainfall is often < 60-100mm (Noy-Meir, 1973). The Tarim River Basin is located in the hinterland of Eurasia, an extremely arid area. This region is highly vulnerable and sensitive to climate change due to the particular features of comprising mountains and a mosaic of inner-mountain basins (Zhou et al., 2020). Populus euphratica (P. euphratica), an ancient, precious, endangered species, is a constructive tree in extremely arid areas (Zhou et al., 2020). P. euphratica is a vital part of the extreme drought ecosystem in this region due to its advantageous characteristics such as cold resistance, heat resistance, drought resistance, wind-sand resistance, and salt-alkali resistance (Li et al., 2019a). In the upper Tarim River basin, since there are less affected by implementing the Ecological Water Conveyance Project than the lower Tarim River, the results that the response of P. euphratica phenology to climate change is more believable. However, research in this eld mainly focuses on the relationship between vegetation and groundwater few if any studies have revealed the dynamics of P. euphratica phenology and its response to climate change in arid areas.
The novelty of this study is that a modi ed method was proposed to explore P. euphratica phenology and its response to climate change using 15-year Global LAnd Surface Satellite (GLASS) leaf area index (LAI) time-series data (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016) in the upper Tarim River basin. A gridded spatial sampling strategy, which can reduce the accuracy error caused by many mixed pixels in satellite remote sensing data, was used to extract the LAI information of P. euphratica from the GLASS LAI time-series data, and the phenological information (i.e., SOS, LOS, EOS) was obtained with the dynamic threshold method from the reconstructed growth time series curve by using the Savitzky-Golay ltering method. The grey relational analysis (GRA) and canonical correlation analysis (CCA) methods were utilized to explore the relationship between phenology and multi-climatic factors. Speci cally, The main objectives of this study were to quantitatively reveal the spatio-temporal characteristics and change trends of P. euphratica phenology, identify the key meteorological factors and key periods that affect the phenology of P. euphratica, and con rm whether the beginning and the end of the P. euphratica growing season in uence each other. The results provide insights into quantitative assessment of spatio-temporal variations of P. euphratica phenology under climate change and provide a scienti c reference for future management of the desert ecosystem.

Study area
Through eld investigation, P. euphratica was found to be mainly distributed along the upper reaches of the Tarim River and the Yerqiang River. Therefore, this investigation mainly focused on the upper Tarim River and the Yeerqiang River Basin, which was divided into three regions (Fig. 1b): Bachu (BC), Alaer (ALE), and Shaya (SY). The study area experiences a typical temperate arid continental climate, with long sunshine duration, strong annual evaporation, and low annual rainfall. The extreme minimum temperature is -25℃, the extreme maximum temperature is 39.4℃, and the annual average temperature

Datasets
The Global LAnd Surface Satellite (GLASS) leaf area index (LAI) has a temporal resolution of eight days and spans the period from 1981-2016. The LAI product from 1982-1999 was developed from AVHRR re ectance with a spatial resolution of 0.05°, and that from 2000-2017 with spatial resolutions of 1 km was derived from MODIS re ectance (Hu et al., 2021). The GLASS LAI data product is a national key project, and its dataset is superior to other similar products in terms of spatial and temporal accuracy (Li and Xiao, 2020). This product feature has signi cant advantages in the continuity of its time series and the integrity of its spatial scope (Xiang et al., 2014). Moreover, the GLASS LAI dataset has passed the veri cation analysis of the measurement data of the common international site, with the results demonstrating that the veri cation accuracy is higher than other similar products (Zhao et al., 2016b;Wang et al., 2017). Since P. euphratica has a small distribution area and low vegetation coverage, this study used the GLASS 1-km resolution data, which were from the National Earth System Science Data Sharing Infrastructure, National Science & Technology Infrastructure of China (http://www.geodata.cn). A total of eight years of data from 2000-2016 were obtained, and all the LAI values were the original values of the GLASS LAI dataset, which were equal to 10 times the actual valued.
This study used climatic data, including sunshine duration (SD), average air temperature (AT), average land surface temperature (ST), effective accumulated temperature (EA), average relative humidity (RH), and precipitation (PR). These data were derived from http://data.cma.cn/ from 2000-2016. Monthly and annual runoff data were derived from the Statistical Bulletin of the Ministry of Water Resources of the People's Republic of China (http://www.mwr.gov.cn/sj/#tjgb). Additionally, in this study, based on Landsat TM/OLI long time series remote sensing data with 30-m resolution from 1990-2015, the distribution data of P. euphratica for 1990, 1995, 2000, 2005, 2010, and 2015 were extracted using the object-oriented classi cation method. These data were then used to select typical sample sites of P. euphratica.

Methods
The technique owchart of this study is shown in Fig. 2, the research steps were as follows. P. euphratica samples were selected and extracted from LAI time series data using the gridded spatial sampling method. The growth time series curve of P. euphratica was then reconstructed using the Savitzky-Golay ltering method. The phenological features of P. euphratica, including the SOS, EOS, and LOS, were extracted using the dynamic threshold method. The Spatiotemporal characteristics of P. euphratica phenology were visualized using the kriging spatial interpolation method and trend analytical method. The grey relational analysis (GRA) method was used to analyze the correlations between climatic factors and the SOS in three periods (spring, SOS preseason, annual), and the EOS in three periods (autumn, EOS preseason, annual). The canonical correlation analysis (CCA) method was used to explore the relationship between phenology (SOS, EOS) and interannual climatic factors (SD, AT, ST, EA, RH, PR).

Extract LAI features
Based on the distribution data of P. euphratica in 1990,1995,2005,2010, and 2015, the regional composition sample data of P. euphratica, which fully met the 1 km 2 area criterion were automatically identi ed by 1 km × 1 km square window. Based on a gridded spatial sampling strategy, the LAI information was extracted from each sample in 2002, 2004, 2006, 2008, 2010, 2012, 2014, and 2016 using the GLASS LAI products.

Filtering and reconstruction of LAI curves
Although the GLASS LAI products have been processed to remove noise, the time series information of the extracted sample points still exhibits the "pseudo" peak phenomenon. Thus, the Savitzky-Golay (S-G) ltering method was utilized to remove the noise and reconstruct the LAI curves of P. euphratica. S-G ltering can be used to the time-series index data of different time scales and is not limited by sensor differences, so the noise in the original data can be removed using the minimum root mean square error (Cai et al., 2017).
Where Y ' j is the resultant LAI value, Y j is the original LAI value, C i is the coe cient for the LAI value of the lter, and N is the number of convoluting integers and is equal to the smoothing window size (2m + 1).

Dynamic threshold method
The dynamic threshold method is widely utilized to extract phenological features of vegetation because it has better spatiotemporal applicability (Hu et al., 2018;Guo et al., 2019). When the dynamic threshold method was proposed to extract phenological features, Jönsson and Eklundh (2002) suggested that the threshold at the SOS and EOS is 0.2, which has been subsequently was used by many scholars in phenological research. This study combined remote sensing data with 2016 ground-based monitoring data and took into account existing research results. The SOS and EOS dates were obtained using various dynamic thresholds ranging from 0.1-0.4 at 0.02 intervals. Then the results were veri ed to use ground-based monitoring data. Based on these tests, we opted to use thresholds of 0.2 and 0.3 to extract the SOS and EOS, respectively. We then used the following procedure, as illustrated by the ow chart in

Grey relational analysis (GRA)
The relationship between many factors falls within a gray area, due to the limitations of people's understanding and the complexity of geographical phenomena, and it is di cult to objectively measure

Canonical correlation analysis (CCA)
Canonical correlation is a statistical analysis method used to study the correlation between two groups of variables and was rst proposed by Hotelling in 1935 (Hotelling, 1935). It is also a dimensionality reduction technology, which is widely used in various elds (Ivanova et al., 2019). Vegetation phenology is not only related to climatic factors, but also the physiological characteristics of vegetation. The change of climatic factors may not lead to one-way phenological changes due to the interaction effects between the SOS in spring and the EOS in autumn (Wu et al., 2016;Shen et al., 2020). It is not rigorous to use simple climatic data to describe these changes. Therefore, the CCA method has been utilized to explore the relationship between phenological features and climatic factors, thereby inferring the interaction between the SOS and EOS. The principle and computation steps are detailed in Hotelling (1935 (Fig. 4d). The slope value was positive, indicating that the SOS showed a trend of postponement with time (Fig. 4e). The SOS in the Bachu (BC) region was delayed by 4.1 days/10a (Fig. 4a), by 5.9 days/10a in the Shaya (SY) region (Fig. 4b), and by 5.9 days/10a in the Alaer (ALE) region (Fig. 4c). The variation range of the SOS in the ALE region was relatively small, while those in the BC and SY regions were larger (Fig. 4e).  (Fig. 5a). The SOS of P. euphratica varied from year to year in the same region and was the earliest in 2004. For the same year in the same region, the SOS dates of P. euphratica in different locations were quite different. In the BC and SY regions, the SOS of P. euphratica exhibited a spatial pattern of early in the west and late in the east (Fig. 5c, d), although there was little difference between the east and west SOS dates in the ALE region (Fig. 5b).

Analysis of EOS
In the study time series, the EOS of P. euphratica showed an advancement trend of about 3.1 days/10a (Fig. 4d). The interannual variation trend of P. euphratica in the SY and ALE regions had a slope < 0, indicating an advancement trend (Fig. 4f), with the EOS advancing at the rate of 5.4 days/10a (Fig. 4b) and 4.6 days/10a (Fig. 4c), respectively. In the BC region, however, the interannual variation trend of the EOS had a slope > 0, indicating a delayed trend (Fig. 4f), at the rate of 0.6 days/10a (Fig. 4a).
The EOS spatial distribution of P. euphratica in the upper Tarim River basin from 2002-2016 is shown in Fig. 6. The end of the P. euphratica growing season in the BC and SY regions became late, while that in the ALE region became earlier (Fig. 6a). The EOS of P. euphratica varied among different years in the same region, and the EOS was the earliest in 2004 (Fig. 6b-d). In the same region, the spatial heterogeneity of the P. euphratica EOS in the same year was higher and the difference was greater. In the ALE and SY regions, the EOS of P. euphratica exhibited a spatial pattern of early in the west and late in the east (Fig. 6a, b), although there was little difference between the east and west EOS dates in the BC region (Fig. 6c).

Analysis of LOS
In the study time series, the LOS of the P. euphratica in the study area displayed a shortening trend, with the LOS shortening 9.2 days/10a (Fig. 7). and the LOS of P. euphratica in the BC, SY, and ALE regions shortening 3.3 days/10a, 12.7 days/10a, and 11.7 days/10a, respectively. The growing season length of P. euphratica differed among the different regions (BC, SY, ALE), with the rank order based on P. euphratica LOS being SY > BC > ALE. The spatial pattern distributions of interannual growth season length of different typical sample sites in the BC and SY regions were different, and the interannual uctuation range was large. The P. euphratica LOS was relatively stable in the ALE region, however. The change trends of the P. euphratica LOS in the study area were similar, and the interannual change trends were the same in the BC, SY, and ALE regions.  (Fig. 8c, 8e). Moreover, the in uence of spring climatic factors in the different regions (BC, ALE, SY) on the SOS of P. euphratica was greater than that of the annual climatic factors or the SOS preseason climatic factors. The EOS of P. euphratica was mainly affected by the meteorological factors (0.69) in autumn (September-November) (Fig. 8b), and less affected by the annual (January-December) meteorological factors (0.52) or the EOS preseason (December-May) meteorological factors (0.67) (Fig. 8d, 8f). Moreover, the in uence of autumn meteorological factors on the EOS of P. euphratica was greater than that of the annual meteorological factors or the EOS preseason meteorological factors in the BC and ALE regions.
Different meteorological factors had different effects on phenology (Fig. 8). On the whole, the sunshine duration (SD) and average air temperature (AT) had the greatest correlation with P. euphratica phenology (SOS, EOS), followed by the average relative humidity (RH), effective accumulated temperature (EA), and average land surface temperature (ST). Precipitation (PR), however, exerted the least effect on P. euphratica phenology. The effects of sunshine duration (0.75), average air temperature (0.76), effective accumulated temperature (0.72), and average relative humidity (0.72) on the SOS of P. euphratica were greater than those of average land surface temperature (0.66) and precipitation (0.30). The effects of sunshine duration (0.86), average air temperature (0.80), effective accumulated temperature (0.75), and average relative humidity (0.74) on the EOS of P. euphratica were greater than those of average land surface temperature (0.68) and precipitation (0.33). Figure 9 shows the CCA patterns between phenology (SOS and EOS) and climatic factors (SD, AT, ST, EA, RH, PR). A total of two pairs of canonical variables were extracted, the rst of which (U 1 , V 1 ) was statistically signi cant (P < 0.05). The results showed that the rst pair of canonical correlation variables revealed a clear correlation between phenology and climatic factors, as indicated by the canonical correlation coe cient of 0.83. The canonical variable U 1 had a signi cant positive correlation with SOS but a signi cant negative correlation with EOS, with typical load coe cients of 0.77 and − 0.98, respectively. The typical variable V 1 had positive correlations with average land surface temperature, average air temperature, and effective accumulated temperature, with the typical load coe cients of 0.71, 0.66, and 0.55, respectively. V 1 had negative correlations with precipitation, sunshine duration, and average relative humidity, however, with typical load coe cients of -0.16, -0.53, -0.30, respectively. These results generally indicated that the SOS had negative correlations with average land surface temperature and average air temperature, and a positive correlation with sunshine duration. The EOS had positive correlations with average land surface temperature and average air temperature, and negative correlations with sunshine duration and average relative humidity.

In uence of interannual climatic change on phenology
The canonical redundancy analysis (Table 1) revealed that the rst pair of canonical variables U 1 could explain 77.62% of the intra-group variation and 18.69% of the other group variation. The typical variable V 1 , however, could explain 27.37% of the intra-group variation and 53.01% of the other group variation.
Therefore, we could infer that the SOS and EOS of P. euphratica affect each other. The phenology (SOS, EOS) of P. euphratica is not only related to meteorological factors but also the physiological characteristics of P. euphratica.  . Based on the annual runoff and rainfall data, this study used the grey correlation method to explore the grey correlation grade between P. euphratica phenology and annual rainfall and runoff. The results demonstrated that the grey correlation grade between runoff and phenology is higher than that between rainfall and phenology (Fig. 11). Precipitation is scarce in arid areas, and P. euphratica growth is not sensitive to precipitation events, which is consistent with previous research results (Xu and Li, 2006;Xu et al., 2011). However, runoff affects the growth of P. euphratica by replenishing the groundwater table, which is then absorbed by the roots of P. euphratica; hence there is a certain lag effect (Liao et al.,2020). Therefore, the grey correlation grade of rainfall and runoff for the phenology of P. euphratica is low, but runoff is still the dominant factor, which is easy to regulate and control. The regulation of runoff and the stable control or reduction of shallow groundwater depth may weaken the impact of global climate change on the phenology of P. euphratica.

Relationships between phenology of P. euphratica and climatic factors
Phenological research of P. euphratica in the Tarim River has been rare. Zhang et al. (2018b) believe that air temperature is the main factor affecting plant phenology in arid and semi-arid areas, and the accumulated temperature should therefore be the condition used to judge P. euphratica phenology. She monitored the phenology of P. euphratica by taking the date the daily average temperature reached ≥ 5 ℃ as the SOS of P. euphratica, and the date of daily average temperature reached ≤ 5 ℃ as the EOS in Hexi Corridor (Zhang et al., 2018b). Her results showed that the SOS of P. euphratica was delayed, while the EOS was advanced, which was consistent with this study's results. However, our study showed that spring sunshine duration (SD) and average air temperature (AT) controlled the SOS and EOS. Moreover, vegetation phenology is affected by many factors such as geographical factors and climatic factors. Using a single accumulated temperature as the criterion for judging the phenology of P. euphratica cannot directly re ect this cycle, and this method is not rigorous.
In this study, the SOS of P. euphratica was delayed by about 5.  2019), who reported that the SOS was advanced and EOS was delayed as the temperature increased. Therefore, it is reasonable that the SOS displays a later trend and EOS shows an earlier trend due to the microclimate of this region is different. In addition, vegetation phenology is regional, which is closely related to species, geographical factors, and climatic factors. The response of vegetation phenology to climate change is not universally identical in different regions (Che et al., 2014;Zhang et al., 2021b). Meanwhile, the spring climatic factors had a great in uence on the SOS of P. euphratica, while the EOS was mainly affected by the autumn climatic factors. As is mentioned above, the response of different vegetation phenology to climate change is not universally identical in different regions. Therefore, it is necessary to explore the variations of P. euphratica phenology and its response to climate change in the extremely arid area.

Advantages and limitations of the study
Through eld investigation, the concentrated distribution of P. euphratica has decreased (Zhao et al., 2016a; Chen et al., 2018). P. euphratica has gradually declined due to climate changes and disturb by human activities, its vegetation coverage is low (Liao et al., 2020; Zhou et al., 2020). the retrieval accuracy of the spatio-temporal patterns of P. euphratica phenology using 1 km x 1 km GLASS LAI remote sensing data is not high. Moreover, considering that there is a wide range of mixed pixels in satellite remote sensing data, the spatial resolution of remote sensing data will affect the accuracy of P. euphratica phenological feature inversion. Therefore, To solve the above two problems, this study utilized a grid spatial sampling strategy to select typical sample areas in order to obtain P. euphratica phenological information. Then, based on P. euphratica spatial distribution data and phenological characteristics data, the Kriging interpolation method was employed to obtain P. euphratica phenological characteristics spatial distribution data in the upper Tarim River basin, thereby reducing the impact of the above two problems. Although there are some errors in the spatial pattern of P. euphratica phenology obtained by the Kriging interpolation method, the study objectives are to explore the changing trend of P. euphratica phenology in arid areas and its response to climate change. Thus, it has little in uence on the spatial distribution pattern and temporal change trend of P. euphratica phenology. This paper has provided an approach that the satellite product data is used to reveal the spatio-temporal variations of the phenology of vegetation, plant community, and a speci c tree species, which the spatial aggregation and coverage are low.
However, due to the limitation of vegetation coverage and aggregation degree of P. euphratica in arid areas, this study needs a longer time-series and higher spatial resolution remote sensing data to reveal the phenological characteristics of P. euphratica. In addition, key periods (spring, SOS preseason, autumn, EOS preseason, annual) that affect the phenology of P. euphratica were divided subjectively, which has a certain subjectivity. Thus, answering these questions await further investigation.

Conclusions
Vegetation phenology is very sensitive to climate change in extremely arid regions, and is thus the best indicator of climate change as well as an essential indicator of ecosystem function. This is of great signi cance to the vegetation ecosystem in arid regions. Using the GLASS LAI time-series data from 2002-2016, the Spatio-temporal characteristics of P. euphratica phenology and their responses to climate change in the upper Tarim River basin were explored. This study revealed that during this period the SOS of P. euphratica was delayed, while the EOS advanced, so the growing season was shortened in the study area. The phenology of P. euphratica was extensively changed, with high spatial heterogeneity. The spring climatic factors had a great in uence on the SOS of P. euphratica, while the EOS was mainly affected by the autumn climatic factors. The ranked order correlation between meteorological factors and the phenology of P. euphratica ranked was SD > AT > EA > RH > ST > PR, implying that it is not easy to characterize phenological change with simple climate data. This di culty is due to the interaction between the SOS in spring and the EOS in autumn, as well as the comprehensive in uence of multiple climatic factors on vegetation phenology change, which have not been adequately considered in current vegetation phenology research and ecosystem models. This study also revealed that the SOS was negatively correlated with ST and AT, but positively correlated with SD, while the relationships between the EOS and meteorological factors were opposite. Due to the scarcity of rainfall in extremely arid regions, vegetation growth primarily depends on groundwater recharge, so the phenology of P. euphratica has a higher correlation with runoff.  do not imply the expression of any opinion whatsoever on the part of Research Square concerning the legal status of any country, territory, city or area or of its authorities, or concerning the delimitation of its frontiers or boundaries. This map has been provided by the authors.     part of Research Square concerning the legal status of any country, territory, city or area or of its authorities, or concerning the delimitation of its frontiers or boundaries. This map has been provided by the authors.