-
PDF
- Split View
-
Views
-
Cite
Cite
Qingling Sun, Jiang Zhu, Siyu Zhu, Baolin Li, Jie Zhu, Xiuzhi Chen, Wenping Yuan, An improved growing season index including the maximum temperature and precipitation to predict foliar phenology of alpine grasslands on the Qinghai–Tibetan Plateau, Journal of Plant Ecology, Volume 18, Issue 2, April 2025, rtaf009, https://doi.org/10.1093/jpe/rtaf009
- Share Icon Share
Abstract
Phenological models are valuable tools for predicting vegetation phenology and investigating the relationships between vegetation dynamics and climate. However, compared to temperate and boreal ecosystems, phenological modeling in alpine regions has received limited attention. In this study, we developed a semi-mechanistic phenological model, the Alpine Growing Season Index (AGSI), which incorporates the differential impacts of daily maximum and minimum air temperatures, as well as the constraints of precipitation and photoperiod, to predict foliar phenology in alpine grasslands on the Qinghai–Tibetan Plateau (QTP). The AGSI model is driven by daily minimum temperature (Tmin), daily maximum temperature (Tmax), precipitation averaged over the previous month (PA), and daily photoperiod (Photo). Based on the AGSI model, we further assessed the impacts of Tmin, Tmax, PA, and Photo on modeling accuracy, and identified the predominant climatic controls over foliar phenology across the entire QTP. Results showed that the AGSI model had higher accuracy than other GSI models. The total root mean square error (RMSE) of predicted leaf onset and offset dates, when evaluated using ground observations, was 12.9 ± 5.7 days, representing a reduction of 10.9%–54.1% compared to other models. The inclusion of Tmax and PA in the AGSI model improved the total modeling accuracy of leaf onset and offset dates by 20.2%. Overall, PA and Tmin showed more critical and extensive constraints on foliar phenology in alpine grasslands. The limiting effect of Tmax was also considerable, particularly during July–November. This study provides a simple and effective tool for predicting foliar phenology in alpine grasslands and evaluating the climatic effects on vegetation phenological development in alpine regions.
摘要
物候模型是预测植被物候及研究植被与气候关系的有效工具。然而,与温带生态系统相比,高寒地区的植被物候模型研究还较少。本研究基于植被生长季指数(GSI),考虑最高气温和最低气温对植被冠层发育的不同影响以及降水和光周期的作用,开发了一个能够模拟和预测整个青藏高原高寒草地叶物候的半经验半机理模型——高寒生长季指数模型(AGSI)。该模型由日最低气温(Tmin)、日最高气温(Tmax)、前一个月平均降水量(PA)和光周期(Photo)4个因子驱动。基于AGSI模型,本研究进一步评估了4个驱动因子对高寒草地叶物候模拟精度的影响,并明确了青藏高原高寒草地返青期和黄枯期的主要气候控制因子。利用地面物候观测数据进行模型评估,结果表明AGSI模型相比其他GSI模型精度更高,AGSI模型模拟的返青期和黄枯期总均方根误差(RMSE)为12.9 ± 5.7天,比其他模型降低了10.9%–54.1%, 包含Tmax和PA的影响使得返青期和黄枯期模拟的总精度提高了20.2%。基于AGSI模型的评估结果表明,4个因子中PA和Tmin对高寒草地叶物候发育的限制作用总体上更加显著和广泛, 而Tmax在7至11月对秋季物候的影响也不容忽视。本研究为预测青藏高原高寒草地叶物候和评估气候对高寒植被物候发育的影响提供了一个简单且有效的工具。
Introduction
Foliar phenology characterizes the start and end of photosynthetic activities of plants and therefore regulates the length of vegetation growing season and terrestrial ecosystem processes, such as carbon and water cycling as well as energy balance between the biosphere and atmosphere (Chambers et al. 2013; Cleland et al. 2007; Richardson et al. 2013). Because the timing of leaf bud burst, unfolding, coloration, and fall is primarily controlled by climatic conditions, plant foliar phenology is highly sensitive to climate change, especially at high latitudes and in high-altitude regions (Park et al. 2020; Suonan et al. 2019; Wu et al. 2022). The Qinghai–Tibetan Plateau (QTP) is the highest plateau in the world and mostly covered by alpine grasslands. It has experienced extraordinary warming at a rate twice the global average during past several decades (Chen et al. 2015a; Yao et al. 2022). Therefore, the phenological shifts of alpine vegetation on the QTP have gained considerable attention over the last few decades (Shen et al. 2015b, 2022).
A large number of studies have explored spatial and temporal shifts of vegetation phenology and linkages of these shifts to climatic factors on the QTP. Previous work is mainly based on remote sensing monitoring (Che et al. 2014; Yu et al. 2010; Zu et al. 2018) and field observation (Sun et al. 2020; Zheng et al. 2016; Zhou et al. 2014). These methods can effectively acquire past and present phenological information, but have distinct limitations in inferring historically missing phenological data, predicting future phenological changes, and investigating the long-term relationship between vegetation and climate (Zhang et al. 2022). Phenological models can make up for these limitations, and there have been a variety of models developed internationally including statistical and process-based models (Piao et al. 2019; Wang et al. 2020a). Statistical models generally use the empirical regression based on environmental factors to predict spring and autumn phenology (Shen et al. 2011; White et al. 1997). They are usually easy to use, but have poor generalization ability. To overcome this deficiency and yield more realistic predictions, many process-based models that explicitly consider plant developmental phases (e.g. endodormancy and ecodormancy) have been developed (Piao et al. 2019). Process-based models are commonly generalizable and have high accuracy, but the prerequisite is clear developmental processes and mechanisms (Zani et al. 2020) . Due to the high requirement for mechanistic processes and simulation accuracy, process-based models are generally calibrated and applied at the species level (Chuine 2000; Delpierre et al. 2009; Fu et al. 2016). There are also some semi-mechanistic models that are neither totally empirical nor completely mechanical, such as the promoter-inhibitor model (Schaber and Badeck 2003), carbon balance model (Arora and Boer 2005), and growing season index model (Jolly et al. 2005). These models are more mechanistic than statistical models, and can be applied at site-to-global scales to obtain both the start and end dates of vegetation growing season. Although a variety of phenological models have been established, most models are developed for woody plants in the temperate and boreal ecosystems (Richardson et al. 2013). Phenological modeling for alpine grasslands has received much less attention.
As for alpine grasslands on the QTP, although few, there are also several phenological models established to estimate spring green-up date or autumn leaf senescence date, including empirical linear models (Shen 2011; Shen et al. 2011), process-based models (Cao et al. 2018; Chen et al. 2015b; Lang et al. 2019), and semi-mechanistic growing season index models (Sun et al. 2018; Wang et al. 2013). Previous efforts on process-based models have been mainly made at plant species level (Chen et al. 2015b; Lang et al. 2019). Cao et al. (2018) have assessed the ability of several state-of-the-art spring phenological models to estimate vegetation green-up dates on the entire QTP and modified existing models by adding environmental constraints. However, no such model is available for predicting vegetation autumn phenology across the entire plateau. The empirical and semi-mechanistic models can be applied at both site and regional scales to estimate foliar phenology of alpine grasslands on the QTP. Nevertheless, these models mainly consider the fixed effect of a single temperature factor, that is usually daily minimum or mean air temperature. It means that the influences of daytime and nighttime temperatures on canopy development are considered the same. However, many previous studies have revealed the asymmetric impacts of daily maximum (or daytime) and minimum (or nighttime) temperatures on both green-up and leaf senescence dates (Fu et al. 2016; Piao et al. 2015; Rossi and Isabel 2017; Wu et al. 2018). Studies on the QTP have also found this phenomenon and shown different effects of increasing daily maximum and minimum temperatures on foliar phenology (Meng et al. 2019; Shen et al. 2016, 2023; Sun et al. 2024; Yang et al. 2017). Previous studies also reported that the influences of temperature on foliar phenology change within the year, particularly the maximum temperature (Liu et al. 2018; Sun et al. 2024; Yu et al. 2010; Zohner et al. 2023). Moreover, there were interactions between air temperature and precipitation on regulating foliar phenology of plants on the QTP (Sun et al. 2023; Suonan et al. 2019). Precipitation plays an important, even a greater, role in determining spring and autumn phenology of alpine grasslands, and also affects the influences of temperature rising on plant phenological shifts (An et al. 2020; Ma et al. 2023; Shen et al. 2015a). Misrepresentation of the complex effects of air temperatures and their interactions with other factors may result in substantial uncertainty when predicting future phenological shifts and ecosystem responses to climate change.
Despite the progresses in phenological models and phenological change mechanism on the QTP, there remains a gap to incorporate the complex effects of daily maximum and minimum temperatures and their interactions with other factors into a model that can be applied to the entire plateau for predicting spring and autumn vegetation phenology. To the best of our knowledge, there has been no such model developed for alpine grasslands on the QTP. There is also no published study that has explicitly evaluated the constraints of multiple influencing factors on foliar development of alpine grasslands and mapped the predominant controls over the entire plateau. Due to high heterogeneity of phenological shifts on the QTP, the model should be flexible to allow primary driving factors to shift in space and time. Growing season index (GSI), a bioclimatic index proposed by Jolly et al. (2005) indicating combined constraints of multiple factors on vegetation canopy development, might meet these requirements. In this study, we tried to establish an improved GSI model incorporating the different effects of daily maximum and minimum temperatures and constraints of preceding precipitation and daylength. Based on the model, leaf onset and offset dates were expected to be predicted more accurately, and the predominant climatic controls on foliar phenology could be identified. To reduce uncertainty brought by data sources and examine the robustness of the improved model, both ground-observed and satellite-derived phenological data were used in this study.
Data and methods
Study area and sites
The Qinghai–Tibetan Plateau (QTP) is the highest region in the world and the largest plateau in Asia. It is mainly located in the western China, and its spatial extent spans from 26°00′ N to 39°47′ N and from 73°26′ E to 104°23′ E, covering a total area of 2.58 × 106 km2 (Fig. 1). Influenced by the Indian Ocean monsoon, westerlies, and terrain, climate on the plateau shows a thermal–moisture gradient from southeast to northwest. The annual mean precipitation decreases from 1764 mm in southeast to 16 mm in northwest during 1981–2011, and the annual mean air temperature decreases from 15.5 °C in southeast to –5.0 °C in northwest (Chen et al. 2015b). Vegetation types mainly consist of grassland, forest, shrub, and cultivated vegetation. The grassland, including steppe and meadow, is dominant and occupies approximately 63.6% of the plateau’s area (Sun et al. 2020). Most plants in grasslands are seasonally deciduous, which have a single growth cycle. The growing season generally ranges from April to October, and the peak occurs during July and August.

Spatial distribution of the 32 study sites (indicated by black circles) and four animal husbandry meteorological sites (indicated by red circles) on the Qinghai–Tibetan Plateau. All the 32 sites are alpine grassland sites, including 11 steppe sites (indicated by outlined circles) and 21 meadow sites (indicated by solid circles).
Meteorological and phenological data at 32 alpine grassland sites were used to calibrate and evaluate the involved GSI models in this study (Supplementary Table S1). The 32 sites include 11 steppe sites and 21 meadow sites, which are all meteorological stations of China Meteorological Administration (CMA) and mainly distributed in the pastoral areas of the QTP. At the 32 sites, satellite phenological data were used to assess the model performance and examine the model robustness. There are four animal husbandry meteorological stations among the 32 sites that have long-term phenological and meteorological observations and plant community survey data, i.e. Haiyan, Henan, Gade, and Qumarleb (Supplementary Table S2). Observed plant species at the four sites are either dominant or widespread species within local grassland communities, with the number of observed plant species ≥ 5 and the length of observed phenological and meteorological data series ≥ 20 (Supplementary Table S2). Thus, foliar phenology on the community scale can be acquired at these four sites through upscaling (see Section Phenological data) for more reliable and comprehensive model evaluations (Supplementary Fig. S1).
Data acquisition and processing
Phenological data
Phenological data used include satellite-derived phenological transition dates and ground-observed plant phenological records. As previously mentioned, long-term and multi-species phenological observation data were available at four sites, i.e. Haiyan, Henan, Gade, and Qumarleb. These data were collected from the annual phenological reports of natural herbaceous plants and growth records of forage from 1983 to 2017. Green-up date (GUD) is defined as leaf onset date, which is identified when 50% of 10 chosen herbs start to display green leaves. Common leaf senescence date (LSD) is defined as leaf offset date, which is identified when leaf coloration in each chosen herb exceeds 2/3. Phenological observations were conducted every two days according to the uniform criterion (China Meteorological Administration 1993). Then, species-level observations were upscaled to obtain the community-level phenological data based on the method of Liang et al. (2011) and Sun et al. (2018) by calculating the weighted means of plant species-specific phenology. The weights of observed species were estimated from their coverage proportions derived from our field plant community surveys (Supplementary Table S2).
Satellite-derived phenological data can provide the start and end dates of grassland growing season at all 32 sites, which have two sources. One is Moderate Resolution Imaging Spectroradiometer (MODIS) Land Cover Dynamics (MCD12Q2) Version 6.1 product, which provides vegetation phenological metrics based on MODIS EVI2 at yearly intervals and 500 m spatial resolution from 2001 to 2020. The other source is Advanced Very High Resolution Radiometer (AVHRR) Vegetation Index and Phenology (VIP) data product (VIPPHEN_EVI2v004), which is provided at 0.05° (ca. 5600 m) spatial resolution and yearly intervals from 1981 to 2014. The two satellite phenological data products were freely downloaded from NASA Land Processes Distributed Active Archive Center (https://disc.gsfc.nasa.gov/), and they have been applied to the QTP (Liu et al. 2024; Zheng et al. 2020). Since alpine grasslands on the QTP have a single growing season in the year, we selected the band “Greenup_1” as start of the season (SOS) and “Dormancy_1” as end of the season (EOS). To acquire effective values from remote sensing data products, geographical locations of the 32 study sites were first corrected using high resolution images on Google Earth. Then, phenological values at these sites were extracted and the same vegetation type was assured in this process.
Meteorological data
We used daily meteorological data observed at the 32 study sites to calibrate and validate all models in this study, including daily maximum air temperature, minimum air temperature, total precipitation, and relative humidity, from 1983 to 2020. These observation data were mainly obtained from China Meteorological Data Service Centre (http://data.cma.cn/), with some additional data from two meteorological bureaus of Qinghai Province and Tibet Autonomous Region in China. Because the meteorological and phenological observations were conducted at almost the same locations in space, these observational data are very useful for modeling and analyzing the phenology–climate relationships (Zheng et al. 2016).
To compare the model performance with previous GSI models, daily photoperiod, VPD, and mean solar radiation were also calculated in this study. Photoperiod was calculated based on site latitude and day of year (Running and Coughlan 1988). VPD was computed using daily maximum air temperature, minimum air temperature, and relative humidity (Tetens 1930). Solar radiation was estimated using MT-CLIM 4.3 (https://www.umt.edu/numerical-terradynamic-simulation-group/project/mt-clim.php), which requires site-related information (i.e. site latitude, elevation, slope, and aspect) as well as daily maximum air temperature, minimum air temperature, and precipitation as input data (Running et al. 1987). Then, the model-derived solar radiation was linearly corrected using ground-observed solar radiation at some sites (Sun et al. 2018).
We also used gridded climatic data to map the predominant climatic controls on foliar phenology of alpine grasslands on the QTP. Daily near-surface air temperature data, including daily maximum and minimum air temperatures, were obtained from National Ecosystem Science Data Center (http://www.nesdc.org.cn). This dataset was released in 2022 and characterized by high spatial resolution and high accuracy (Fang et al. 2022). Daily gridded precipitation data were collected from National Tibetan Plateau Data Center (http://data.tpdc.ac.cn). It is a new gridded precipitation dataset for Chinese mainland constructed based on gauge observations, and has been demonstrated to be a high-quality dataset (Han et al. 2023). The spatial resolution of these climatic data was 0.1° and data during 1979–2018 were used to calculate the multi-year average temperature and precipitation indices in the AGSI model.
Vegetation index data
The normalized difference vegetation index (NDVI) is commonly used to represent canopy greenness dynamics. We used two NDVI products, MODIS NDVI and SPOT-VGT NDVI, to explore the relationship between AGSI and NDVI. MODIS NDVI used is MOD13A1 version 6.1, derived from NASA Distributed Active Archive Center (https://ladsweb.modaps.eosdis.nasa.gov/search/). It provides NDVI values at 500 m spatial resolution and the temporal resolution is 16 days. SPOT-VGT NDVI was derived from VITO VEGETATION center (http://free.vgt.vito.be), having a spatial resolution of 1 km and temporal resolution of 10 days. For comparison at each site, NDVI values closest to each site were first extracted based on site latitude and longitude. Then, we calculated multi-year average NDVI and daily AGSI values during the same time periods using MODIS NDVI data from 2000 to 2017 and SPOT NDVI data from 2000 to 2012. At last, we identified the maximum AGSI in the corresponding 16-day or 10-day period of NDVI data, and compared the maximum AGSI with NDVI values using a standard Pearson’s correlation.
Vegetation type data
The vegetation distribution and grassland type information were obtained from a digitized 1:1 000 000 vegetation map of China (Chinese Academy of Sciences, 2001). The grassland in the vegetation map has two categories, steppe and meadow (Fig. 1). We further extracted and aggregated the areas of alpine meadows and alpine steppes according to the subcategories of meadow and steppe on the QTP, and then the extent of alpine grasslands could be obtained. During our study period, we assumed that there was no change in the grassland distribution and type.
Model development
We tried to develop a simple and effective prognostic phenology model considering the effects and interactions of daily maximum temperature, daily minimum temperature, precipitation, and daylength, for alpine grasslands on the QTP under the framework of GSI. GSI is originally derived from three environmental factors, i.e. daily minimum air temperature, mean daily vapor pressure deficit (VPD), and daily photoperiod (Jolly et al. 2005). The three factors are further converted to three indices through linear functions ranging from 0 to 1, indicating temperature, humidity, and light constraints on canopy development, respectively. By multiplying these three indices, daily iGSI is formed, representing a concurrent control of temperature, humidity, and light on foliar phenology. The last step is to calculate the 21-day moving average of iGSI and use a threshold to identify leaf onset and offset dates. The moving average serves to buffer extreme events from prematurely triggering canopy changes.
The effectiveness and flexibility of the GSI model have been recognized worldwide (Hidy et al. 2012; Orlandi et al. 2013; Stöckli et al. 2008). The GSI model has also been applied to the QTP and improved for alpine meadows (Sun et al. 2018; Wang et al. 2013). However, existing GSI models only used daily minimum temperature to indicate low temperature constraint, ignoring the important and different effect of daily maximum temperature. Here, we added an additional index into GSI to indicate the impact of daily maximum temperature on foliar development of alpine grasslands on the QTP. In addition, we adopted a precipitation index, mean precipitation over the previous month, to replace the VPD in the original GSI. This is because precipitation plays a critical role in regulating both spring and autumn foliar phenology (Ma et al. 2023; Shen et al. 2015a), and it shows better performance in phenological modeling (Sun et al. 2018). The improved GSI model was driven by four factors: daily minimum air temperature (Tmin), daily maximum air temperature (Tmax), precipitation averaged over the previous month (PA), and daily photoperiod (Photo).
Similarly, daily indices that varied between 0 (indicating totally constrained) and 1 (indicating totally unconstrained) were first calculated. There were totally four indices, hereafter termed as iTmin, iTmax, iPA, and iPhoto, corresponding to the Tmin, Tmax, PA, and Photo, respectively. Then, an improved GSI for alpine grasslands, referred to as the AGSI, was computed as the 21-day running average of the daily index iAGSI that was the product of iTmin, iTmax, iPA, and iPhoto (Eq. 1). A higher value of AGSI indicates weaker constraints on canopy development and vice versa. The last step was to identify leaf onset and offset dates using an AGSI threshold (Fig. 2).

Sketch map showing the identification of leaf onset date (indicated by A) and leaf offset date (indicated by B) based on the AGSI curve and an AGSI threshold (AGSIthr). Leaf onset and offset dates were identified as the time when the smoothed AGSI curve (i.e. 21-day moving average of iAGSI) was first and last above the AGSIthr, respectively.
iTmin
The main manifestation of the effect of daily minimum temperature on vegetation foliar development is low temperature constraint, which is closely related to the frost and freezing damage as well as vegetation respiratory capacity (Meng et al. 2019; Shen et al. 2016). We assumed that very low daily minimum temperature was lethal, which could fully restrict plant growth and development. Daily minimum temperature higher than a threshold was assumed to have no constraint on foliar phenology. Calculation formula of iTmin was kept the same as that in Jolly et al. (2005) (Eq. 2),
where Tmin is daily minimum air temperature (°C). Tmmin and Tmmax are the lower and upper minimum temperature thresholds, respectively.
iTmax
The effect of daily maximum temperature on foliar development is complex. We assumed that daily maximum temperature lower than a threshold could fully restrict foliar development due to very small photosynthesis and evapotranspiration rates under extremely low daytime temperatures (Lange et al. 1981). On the other hand, extremely high daily maximum temperatures could also restrain photosynthesis and cause water stress, which may lead to more allocation to the underground parts and be unfavorable for foliar development (Zhu et al. 2023). So, daily maximum temperature higher than another threshold could also totally constrain canopy development. The optimum air temperatures ranging from 15 °C to 20 °C were defined for alpine grasslands on the QTP (Cui et al. 2012; Huang et al. 2019; Zhao 2009), and we assumed that daily maximum temperatures within this range had no limiting effect on foliar development. Therefore, we used a piecewise function to calculate iTmax as follows (Eq. 3),
where Tmax is daily maximum air temperature (°C). Tamin and Tamax are the lower and upper maximum temperature thresholds, respectively.
iPA
Similar to Sun et al. (2018), we used accumulated precipitation to represent the main source of soil water and average precipitation over the previous month as a surrogate for mean soil moisture, considering that daily precipitation is discrete and soil moisture is relatively stable. The period of one month was determined because previous studies had showed that precipitation during one previous month was most relevant to the timing of green-up and leaf senescence on the QTP (Ma et al. 2023; Piao et al. 2006; Sun et al. 2024). We also conducted an experiment to find the optimum time period for precipitation accumulation leading to the lowest modeling error with a step of 5 days, and this period was also around 30 days (Fig. S2). iPA was calculated as follows,
where PA is precipitation averaged over one previous month (mm), and P is the daily precipitation (mm). The d is day of year, and if d is lower than 30, precipitation data from the last year will be used to calculate PA. PAmin and PAmax are the lower and upper thresholds of PA, respectively.
iPhoto
Photoperiod is a reliable annual climatic cue because it does not vary from year to year at a given location. The effect of photoperiod is that it provides the outer envelope within which other climatic factors may dictate foliar development (Jolly et al. 2005). Studies have shown that photoperiod is important to both leaf flush and senescence, and interacts with temperature to limit foliar phenology (Delpierre et al. 2009; Fu et al. 2019). iPhoto was calculated as Jolly et al. (2005) (Eq. 6),
where Photo is daily photoperiod (h). Photomin and Photomax are the lower and upper photoperiod thresholds, respectively.
Model calibration and validation
Due to a relatively small date size, we used 10 × 3 fold cross-validation to evaluate model performance, which means that three-fold cross-validation was implemented 10 times. Each time the phenological and meteorological datasets were divided into three parts, and two-thirds parts were used to calibrate the models (calibration dataset) and the remaining part was used to validate the models (validation dataset). Because we used 10 different random seeds to divide the whole dataset, 30 different calibration and validation datasets could be obtained by using the 10 × 3 fold cross-validation method. Accordingly, model calibration was conducted 30 times using 30 different calibration datasets, which generated 30 sets of optimized model parameters. Model validation was also conducted 30 times using the 30 different validation datasets and optimized model parameters. The average evaluation indices derived from the 30 times validation were used to quantitatively assess model performance.
The dual annealing optimization method was used to implement model calibration because this algorithm demonstrated the best performance among several non-linear optimization algorithms (Zhdanov et al. 2023). The dual annealing, as implemented by the scipy.optimiza.dual_annealing function of the scipy toolbox, combines the classical simulated annealing with fast simulated annealing algorithms augmented by a local search on accepted locations (Xiang et al. 1997). Parameters needed to be optimized included Tmmin, Tmmax, Tamin, Tamax, PAmin, PAmax, Photomin, Photomax, and an AGSI threshold (AGSIthr). The ranges and initial values of these parameters were determined based on previous studies (Jolly et al. 2005; Stöckli et al. 2008; Sun et al. 2018) and local climatic observations (Supplementary Table S3). The default values for the maximum number of global iterations (1000) and the limit for the number of objective function calls (107) were used. The only constraint was that optimized values indicating lower thresholds should be lower than those indicating upper thresholds. The energy function applied in the optimization was the lowest total root mean square errors (RMSE) between the predicted and observed leaf onset and offset dates (Eq. 7),
where E is the total RMSE and n is the number of years. and indicate ground-observed/satellite-derived leaf onset and offset dates in the ith year, respectively. and indicate model-predicted leaf onset and offset dates in the ith year, respectively.
We also compared performance of the AGSI model with four existing GSI models that had been applied to grasslands worldwide. These models include the original GSI model proposed by Jolly et al. (2005) (termed as GSI0), the HSGSI model proposed by Hidy et al. (2012) for temperate grasslands, a revised GSI model that had been applied to an alpine meadow (Wang et al. 2013; referred to as WGSI), and the AMGSI model proposed by Sun et al. (2018) for alpine meadows on the QTP. None of the four models has taken into account different impacts of daily maximum and minimum temperatures, and only the AMGSI model has considered the effect of precipitation. All the models have used daily minimum air temperature to reflect the low temperature constraint and daily photoperiod to represent the sunshine condition, except the WGSI that replaces photoperiod with global solar radiation. More detailed descriptions on these models are shown in Table 1.
Detailed descriptions on the GSI models used in this study. The performance of the AGSI model was compared with that of other four existing models
Model . | Key model equation . | Number of parameters . | Reference . |
---|---|---|---|
GSI0 | iGSI0 = iTmin × iVPD × iPhoto | 7 | Jolly et al. 2005 |
HSGSI | iHSGSI = iTHS10 × | 9 | Hidy et al. 2012 |
WGSI | iWGSI = iTmin × iVPD × iRad | 7 | Wang et al. 2013 |
AMGSI | iAMGSI = iTmin × iPA × iPhoto × iRad × iSnow | 9 | Sun et al. 2018 |
AGSI | iAGSI = iTmin × iTmax × iPA × iPhoto | 9 | This study |
Model . | Key model equation . | Number of parameters . | Reference . |
---|---|---|---|
GSI0 | iGSI0 = iTmin × iVPD × iPhoto | 7 | Jolly et al. 2005 |
HSGSI | iHSGSI = iTHS10 × | 9 | Hidy et al. 2012 |
WGSI | iWGSI = iTmin × iVPD × iRad | 7 | Wang et al. 2013 |
AMGSI | iAMGSI = iTmin × iPA × iPhoto × iRad × iSnow | 9 | Sun et al. 2018 |
AGSI | iAGSI = iTmin × iTmax × iPA × iPhoto | 9 | This study |
Note: Tmin and Tmax are daily minimum and maximum air temperatures, respectively; VPD is mean daily vapor pressure deficit; Photo is daily photoperiod; HS10 is 10-day heatsum calculated as the sum of mean daily temperature−5 °C for each 10 day-long period; Rad is daylight mean global solar radiation; PA is the precipitation averaged over the previous month; Snow is a boolean index indicating whether there is snowfall based on mean daily temperature and precipitation data.
Detailed descriptions on the GSI models used in this study. The performance of the AGSI model was compared with that of other four existing models
Model . | Key model equation . | Number of parameters . | Reference . |
---|---|---|---|
GSI0 | iGSI0 = iTmin × iVPD × iPhoto | 7 | Jolly et al. 2005 |
HSGSI | iHSGSI = iTHS10 × | 9 | Hidy et al. 2012 |
WGSI | iWGSI = iTmin × iVPD × iRad | 7 | Wang et al. 2013 |
AMGSI | iAMGSI = iTmin × iPA × iPhoto × iRad × iSnow | 9 | Sun et al. 2018 |
AGSI | iAGSI = iTmin × iTmax × iPA × iPhoto | 9 | This study |
Model . | Key model equation . | Number of parameters . | Reference . |
---|---|---|---|
GSI0 | iGSI0 = iTmin × iVPD × iPhoto | 7 | Jolly et al. 2005 |
HSGSI | iHSGSI = iTHS10 × | 9 | Hidy et al. 2012 |
WGSI | iWGSI = iTmin × iVPD × iRad | 7 | Wang et al. 2013 |
AMGSI | iAMGSI = iTmin × iPA × iPhoto × iRad × iSnow | 9 | Sun et al. 2018 |
AGSI | iAGSI = iTmin × iTmax × iPA × iPhoto | 9 | This study |
Note: Tmin and Tmax are daily minimum and maximum air temperatures, respectively; VPD is mean daily vapor pressure deficit; Photo is daily photoperiod; HS10 is 10-day heatsum calculated as the sum of mean daily temperature−5 °C for each 10 day-long period; Rad is daylight mean global solar radiation; PA is the precipitation averaged over the previous month; Snow is a boolean index indicating whether there is snowfall based on mean daily temperature and precipitation data.
Model performance was assessed mainly using two indices, i.e. root mean square error (RMSE) and mean error (ME). The RMSE is used as an accuracy metric (Lang et al. 2019; Tao et al. 2021), which measures the absolute differences between values predicted by a model and values that are actually observed (Eq. 8). The ME is a bias metric (Descals et al. 2020), which is defined as the mean difference between model predictions and observations (Eq. 9). A positive ME indicates that the model-predicted leaf onset/offset dates were later than the observations, and vice versa. The best model was the one having the lowest RMSE and the ME closest to zero.
Here, OBSi and MODi represent the observed and modeled leaf onset/offset dates in the ith year, respectively. n is the number of observations.
Contribution assessment of individual factors to modeling accuracy
To quantitatively assess the contributions of individual factors in the AGSI model to predictions of foliar phenology, we removed one factor and the corresponding index (indicating the limiting effect of this factor) each time from the AGSI model to obtain new formulations (Supplementary Table S4). At each site, the AGSI model was calibrated and validated using dual annealing algorithm and 10 × 3 fold cross-validation. A total of 30 sets of parameter values were obtained, and among them, one set of parameter values that minimized the total RMSE of predicted leaf onset and offset dates was selected as the optimum parameter values. Based on these optimum parameter values, all daily indices in the AGSI model were calculated. The contribution of a factor to modeling accuracy, which reflects the influence of the factor on phenological prediction, was evaluated by the increases in the predicted RMSEs after removing this factor and corresponding index from the AGSI formulation (Eq. 1; Supplementary Table S4). In the process of removal, only one factor was eliminated each time, and the remaining factors and their indices were kept unchanged.
Investigation of predominant climatic controls on foliar phenology
Based on the AGSI model, we further identified and explicitly mapped the most predominant climatic controls on foliar phenology over the entire alpine grasslands on the QTP. Because photoperiod has no inter-annual variation, and strictly speaking, it is not a climatic factor, its limiting effects were not compared with those of temperatures and precipitation. Daily gridded air temperature and precipitation data from 1979 to 2018 were used to calculate the point-wise daily indices of iTmin, iTmax, and iPA (Eqs. 2–5) using model parameter values calibrated from ground observations (ground-based) and MODIS phenological data (satellite-based). Two different sets of parameter values were obtained for alpine meadows and steppes through model calibration (Supplementary Table S5). Considering that all daily indices varied between 0 (totally constrained) and 1 (totally unconstrained), we first subtracted the multi-year average iTmin, iTmax, and iPA from 1.0 to represent daily limiting effects of Tmin, Tmax, and PA, respectively. Then, the daily constraints were summed over January–May (time period relevant to leaf onset dates) and July–November (relevant to leaf offset dates) to indicate the total constraints on leaf onset and offset, respectively (Eqs. 10–12). The most predominant controls on leaf onset and offset dates were identified as the factors that showed the maximum limiting effects during January–May and July–November, respectively (Eq. 13).
Here, t1 and t2 are the start and end days of year to accumulate the limiting effects of climatic factors. t1 and t2 were set to 1 and 153, respectively, when calculating the climatic limiting effects on leaf onset dates. Otherwise, their values were set to 182 and 334 for calculating the constraints on leaf offset dates. LETmin, LETmax, and LEPA, indicate the total limiting effects of Tmin, Tmax and PA, respectively.
Results
Evaluation of model performance
Based on ground-observed climatic and phenological data, individual daily index values and the AGSI curves could be obtained at each site (Supplementary Fig. S3). The multi-year average AGSI and MODIS/SPOT NDVI at the same sites showed very significant positive correlations (P < 0.001), which indicates that model-predicted AGSI well represented canopy greenness dynamics at all sites regardless of the dominant climatic controls at that site (Fig. 3). Correlations between AGSI and SPOT NDVI were better than those between AGSI and MODIS NDVI because of a higher temporal resolution of SPOT NDVI data. The AGSI and NDVI curves showed similar seasonal patterns on the whole, but there were also inevitable inconsistencies between the AGSI and NDVI curves in some cases. For example, the average AGSI and NDVI curves at Qumarleb and Henan began to show clear differences since the peaks of AGSI curves. In the former half of the growing season, the AGSI curves fit NDVI curves better, indicating that the AGSI model is conducive to better predictions of leaf onset dates at these sites.

Comparison of the multi-year average AGSI and NDVI at four sites. The NDVI curves were acquired using (a, c, e, g) MODIS NDVI from 2000 to 2017 and (b, d, f, h) SPOT NDVI from 2000 to 2012, which have a temporal resolution of 16 days (n = 23) and 10 days (n = 36), respectively. AGSI values were the maximum values in the corresponding 16-day or 10-day period of NDVI data. The r and P indicate the Pearson’s correlation coefficient and significance level, respectively.
There is an overall good consistence between the simulated leaf onset and offset dates with observed green-up and leaf senescence dates (Fig. 4). The predicted leaf onset and offset dates were highly correlated with ground observations at the four sites, with R2 of 0.75 and 0.69, respectively. The regression line between simulated leaf onset date and observed green-up date was very close to the 1:1 line and the scatters were distributed uniformly on both sides of the 1:1 line. Slope of the regression line between simulated leaf offset date and observed leaf senescence date was smaller, because of the relatively large differences at Qumarleb and Henan. It means that the AGSI model had overall better performance in simulating the interannual variability of leaf onset dates. Comparison among the four sites showed better consistency between the simulated and observed dates at Gade and Haiyan sites (Supplementary Fig. S4).

Comparison of the (a) simulated leaf onset dates and observed green-up dates, and (b) simulated leaf offset dates and observed leaf senescence dates at the four animal husbandry sites on the QTP, including Qumarleb, Gade, Henan, and Haiyan. Each dot represents the simulated and observed dates in one year at a site, which is Julian day of the year. The dashed line is 1:1 line, and the solid line is regression line. The confidence intervals are shaded in gray.
Compared with other GSI models, the AGSI model showed higher accuracy in predicting both leaf onset and offset dates when assessed using ground observations (Fig. 5a). On average, the RMSEs of leaf onset and offset dates predicted by the AGSI model were 6.3 d and 6.6 d, respectively, decreasing those of existing GSI models by 15.8%–43.9% and 6.0%–60.9%, respectively. The total RMSE of predicted leaf onset and offset dates based on the AGSI model was 12.9 days. The lowest total RMSE calculated from other models was 14.6 days of the AMGSI, and the highest total RMSE was 28.3 days of the WGSI. The AGSI model decreased total RMSE by 10.9%–54.1% compared with other GSI models. Compared with the GSI0, RMSEs of leaf onset and offset dates predicted by the AGSI model decreased by 22.8% and 17.5%, respectively. The total modeling accuracy of leaf onset and offset dates was improved by approximately 20.2% after considering the effects of Tmax and PA. The MEs of predicted leaf onset and offset dates from the AGSI model were 0.37 day and 0.27 day, respectively, which indicates that model-predicted leaf onset and offset dates were slightly later than the observed dates (Fig. 5b). Although the MEs from the GSI0 and AMGSI models were also small, their values were negative, indicating that model-predicted leaf onset and offset dates were earlier than the observed dates. The MEs of the WGSI and HSGSI models were clearly higher, suggesting distinct biases in the predicted leaf onset and offset dates.

Comparison of model performance evaluated by the (a) RMSE and (b) ME using ground observations at the four animal husbandry meteorological sites. AGSI was the model proposed in this study. GSI0, HSGSI, WGSI, and AMGSI were four existing GSI models that had been applied to grasslands. Averages show the mean results of the four sites.
We noticed that model performance varied widely between different sites (Fig. 5). The total RMSEs of the AGSI model were 19.5 days, 8.5 days, 15.9 days, and 7.9 days at Qumarleb, Gade, Henan, and Haiyan, respectively. The total MEs ranged from –0.6 days at Gade to 1.6 days at Qumarleb. The AGSI model showed higher accuracy in predicting leaf onset dates at Qumarleb and Henan, while higher accuracy in modeling leaf offset dates was found at Gade and Haiyan. The MEs of predicted leaf onset and offset dates at Qumarleb and Haiyan were all positive, while those at Gade and Henan had negative values. The RMSEs and MEs of other models also varied distinctly between sites. Despite the distinctions in model performance at different sites, the AGSI model always showed the lowest RMSE among all GSI models, and the MEs were all close to zero. This indicates that the AGSI model is robust in space.
Evaluation of model performance using different satellite phenological data showed similar results (Fig. 6). The AGSI model showed the lowest total RMSEs of predicted leaf onset and offset dates, which were, on average, 24.5 days and 29.8 days when assessed using MODIS and AVHRR phenological data at all alpine grassland sites, respectively. Compared with other models, the AGSI model decreased the total RMSEs of predicted leaf onset and offset dates by 1.5–14.5 days (assessed using MODIS data) and 2.0–13.5 days (assessed using AVHRR data), respectively. Through model comparison between the AGSI and GSI0, we found that the total modeling accuracy of leaf onset and offset dates was improved by 22.9% (assessed using MODIS data) and 17.2% (assessed using AVHRR data), respectively, after considering Tmax and PA in the AGSI. All the MEs of modeled leaf onset and offset dates were negative regardless of models and satellite sources, which indicates that model-predicted leaf onset and offset dates were earlier than satellite-derived SOS and EOS. The total MEs from the AGSI model were –10.1 days and –20.7 days when assessed using MODIS and AVHRR data, respectively, and the biases were smaller than those of other GSI models. The model evaluation results using two different types of satellite phenological data also demonstrated the robustness of the AGSI model.

Comparison of model performance indicated by the RMSE and ME using satellite phenological data at (a, d, g) alpine grassland, (b, e, h) alpine meadow, and (c, f, i) alpine steppe sites. The dots and diamonds indicate averages acquired using MODIS and AVHRR phenological data, respectively. The horizontal and vertical error bars represent one standard deviation in the RMSE and ME, respectively.
Based on satellite phenological data, we also found that model performance varied between different sites, and the AGSI model had overall higher accuracy at alpine meadow sites (Fig. 6). The total RMSEs of the AGSI model calibrated and validated using MODIS phenological data were on average 23.2 days and 27.5 days at alpine meadow and steppe sites, respectively. The mean total RMSEs acquired from AVHRR VIPPHEN data were 29.3 days and 31.9 days at alpine meadow and steppe sites, respectively. Despite higher modeling accuracy at alpine meadow sites, the biases in model predictions at alpine meadow sites were not necessarily smaller than those at alpine steppe sites. The total MEs of the AGSI model assessed using MODIS MCD12Q2 data were, on average, –9.7 days and –10.4 days at alpine meadow and steppe sites, respectively. However, when evaluated using AVHRR VIPPHEN data, the mean total MEs at alpine meadow and steppe sites were –20.6 days and –15.1 days, respectively.
Contributions of different factors to modeling accuracy
Based on ground-observed data at four sites, we found that all four factors in the AGSI had important impacts on predictions of foliar phenology (Supplementary Table S6). Excluding Tmax, Tmin, PA, and Photo from the AGSI model, the total RMSEs of predicted leaf onset and offset dates increased by 5.5 days (53.4%), 7.1 days (68.9%), 10.1 days (98.1%), and 27.4 days (266.0%), respectively. Photoperiod showed an overall highest contribution mainly due to its critical role in predicting leaf offset dates. The RMSEs of predicted leaf onset and offset dates increased by 4.3 days and 23.1 days, respectively, after removing iPhoto from the AGSI model. Precipitation was secondary to photoperiod, whereas it played a more important role in modeling leaf onset dates. Excluding iPA from the AGSI, RMSEs of predicted leaf onset and offset dates increased by 7.7 days and 2.4 days, respectively. Tmin also showed a greater impact on leaf onset dates, and without iTmin, RMSEs of predicted leaf onset and offset dates increased by 5.0 days and 2.1 days, respectively. Tmax also had considerable effects because RMSEs of predicted leaf onset and offset dates increased by 2.7 days (50.9%) and 2.8 days (56.0%), respectively, after removing iTmax from the AGSI model.
Contribution assessment results varied between different sites (Fig. 7). At Qumarleb, precipitation showed a more important influence than daily maximum and minimum temperatures on modeling foliar phenology. Although their impacts on predicting leaf offset dates were similar, precipitation played a more critical role in modeling leaf onset dates. At Gade, the contributions of daily maximum temperature outweighed those of daily minimum temperature and precipitation, particularly on predictions of leaf onset dates. At Henan, the overall impact of precipitation on predicting foliar phenology was almost equal to that of minimum temperature, and larger than that of maximum temperature. At Haiyan, an alpine steppe site, precipitation had a more distinct contribution to accurate predictions of leaf onset dates than temperatures. Despite large differences, photoperiod was always the most critical factor determining model accuracy at all sites.

Impacts of individual factors in the AGSI model on prediction accuracy at different sites, which were assessed using ground observations. AGSI–iTmax, AGSI–iTmin, AGSI–iPA, and AGSI–iPhoto represent the revised AGSI formulations excluding iTmax, iTmin, iPA, and iPhoto from the AGSI, respectively.
Evaluations using two sources of satellite phenological data showed broadly similar results (Fig. 8). Overall, precipitation had the most prominent influences on modeling accuracy, especially on predictions of leaf onset dates. Impacts of daily minimum and maximum temperatures were also considerable, yet not greater than precipitation. Interestingly, photoperiod, the most critical factor affecting ground-observed leaf offset dates, showed much smaller effects when evaluated using satellite phenological data. All the factors appeared to well regulate leaf onset dates, but showed not much effects on leaf offset dates. The mean RMSE of predicted leaf offset dates increased slightly, even decreased, after removing a factor from the AGSI model. As regards leaf onset dates, precipitation showed the most predominant effects. Once excluding precipitation from the AGSI model, mean RMSE of modeled leaf onset dates increased by 45.3 days and 30.9 days when evaluated using MODIS and AVHRR phenological data, respectively. Compared–with the maximum temperature, minimum temperature had a greater contribution to predictions of leaf onset dates, but their impacts on leaf offset dates were very close.

Impacts of Tmax, Tmin, PA, and Photo on accuracy of the AGSI model, which were assessed using satellite phenological data at (a–b) all alpine grassland sites, (c–d) alpine meadow sites, and (e–f) alpine steppe sites, respectively. The box plots give the 25th (Q1) and 75th (Q3) quantile, and middle of each box is the 50th quantile. The black dots inside boxes indicate mean values. The lower and upper whiskers are Q3 + 1.5IQR and Q1–1.5IQR, respectively, with IQR = Q3–Q1.
At both alpine meadow and steppe sites, the four factors in the AGSI model could effectively regulate leaf onset dates because the RMSEs of modeled leaf onset dates increased significantly after excluding a factor from the AGSI (Fig. 8). However, they had not much impacts on predictions of leaf offset dates, indicating that all four factors had limited influences on leaf offset dates at both alpine meadow and steppe sites. Although the impacts of individual factors were overall similar, there were still some differences between alpine meadow and steppe sites. The impacts of minimum temperature on modeling leaf onset dates at alpine meadow sites seemed to be larger than those at alpine steppe sites. Instead, precipitation and maximum temperature had greater effects at alpine steppe sites. The overall effect of photoperiod was minimal, but it had greater impacts at alpine meadow sites by comparison.
Predominant climatic controls on foliar phenology
Using parameter values calibrated from ground observation data, the AGSI model showed that precipitation rather than temperatures was the most predominant climatic control on foliar development of alpine grasslands on the QTP (Fig. 9a–f). From January to May, there were 61.0% and 37.6% of QTP alpine grasslands in which foliar phenology (i.e. leaf onset dates) was most limited by precipitation and the minimum temperature, respectively. Only 1.4% of alpine grasslands were most restricted by daily maximum temperature. The limiting effects of precipitation were greater than other factors in approximately 62.5% of alpine meadows and 59.6% of alpine steppes. Constraints of minimum temperature were stronger in 35.0% of alpine meadows and 39.8% of alpine steppes. The predominant controls of precipitation and minimum temperature were widely distributed across the plateau, but the limiting effects of minimum temperature were more prominent in the northern and southwestern regions of alpine grasslands where air temperature was relatively low.

Spatial distribution of the most predominant climatic controls on leaf onset and offset dates of alpine grasslands on the QTP, identified based on the AGSI model using (a–f) ground-based and (g–l) satellite-based parameter values, respectively. The most predominant climatic controls on leaf onset and offset dates were identified using daily indices from January to May and July to November, respectively.
From July to November, precipitation and minimum temperature showed a more important influence on foliar phenology (i.e. leaf offset dates) in 56.3% and 29.1% of alpine grasslands, respectively (Fig. 9d–f). The influences of precipitation were more prominent than other factors in approximately 48.9% of alpine meadows and 62.9% of alpine steppes. Constraints of daily minimum temperature were greatest in 29.7% of alpine meadows and 28.6% of alpine steppes. The limiting effects of daily maximum temperature were also considerable, and about 14.6% of QTP alpine grasslands were restricted predominantly by the maximum temperature. The restrictions of maximum temperature were more extensive in alpine meadows because there were 21.3% of alpine meadows limited primarily by maximum temperature during July–November. The constraints of temperatures were mainly concentrated in the northeastern plateau (relatively moist), and the limiting effects of precipitation were more prominent in the central and western parts (relatively arid).
Using satellite-based parameter values, modeled predominant climatic controls were also mainly precipitation and minimum temperature (Fig. 9g–l). There were 53.6% and 46.1% of QTP alpine grasslands in which foliar development from January to May was most restricted by precipitation and minimum temperature, respectively. Precipitation showed the most extensive importance in predicting leaf onset dates of alpine steppes because there were 68.5% of alpine steppes limited predominantly by precipitation during January–May, which were concentrated in the western plateau (relatively arid). As for alpine meadows, minimum temperature showed the greatest limiting effects on leaf onset dates because 62.7% of alpine meadows were restricted predominantly by the minimum temperature during January–May, which were mainly distributed in the central, southeastern, and southwestern regions (relatively moist and cold). The effects of maximum temperature on leaf onset dates seemed to be slight compared to those of precipitation and minimum temperature.
Satellite-based model investigation also showed that precipitation was the most extensive limiting factor of leaf offset dates, followed by the minimum temperature (Fig. 9j–l). Overall, there were 56.7% and 38.3% of QTP alpine grasslands in which foliar development from July to November was most limited by precipitation and minimum temperature, respectively. About 5.0% of alpine grasslands were dominated by limitations of maximum temperature, which were mainly distributed in the eastern and northeastern plateau. From July to November, there were 47.6% of alpine meadows and 64.8% of alpine steppes limited predominantly by precipitation. There were 46.2% of alpine meadows and 31.4% of alpine steppes restricted most by minimum temperature. Although the constraints of maximum temperature during July–November were smaller than those of precipitation and minimum temperature, they were much stronger than the effects during January–May. We also found that the effects of precipitation in alpine steppes were greater than those in alpine meadows, whereas the constraints of air temperatures were just opposite.
Discussion
Modeling differences between ground and satellite phenological data
Performance of the AGSI model assessed using ground-observed phenological data was substantially better than that evaluated using satellite-derived phenological data. Comparison between two satellite phenological data showed that the model accuracy was higher when using high-resolution MODIS data than using coarser AVHRR data. This suggests that phenological data used in modeling could affect accuracy, and higher accuracy could be expected with smaller data uncertainty and higher data resolution. As regards the impacts of individual factors in the AGSI, evaluation results based on satellite phenological data were significantly different from those using ground phenological data. There were two main differences. One is that satellite-based evaluations showed significantly weaker effects of photoperiod and daily maximum temperature on predictions of foliar phenology compared with ground-based results. The other difference is that, although both ground-based and satellite-based evaluations found greater impacts of the four factors in the AGSI model on predictions of leaf onset dates, satellite-based evaluations showed much smaller effects of the four factors on leaf offset dates by comparison.
Above discrepancies can be attributed to the contrasting properties and changes of ground-observed and satellite-derived phenological data. Ground-observed phenology is the traditional plant phenological transition dates, which indicates the development stage of plant individuals at local scales (Sun et al. 2020). Satellite-derived phenology is land surface phenology, indicating changes in vegetation greenness at landscape or regional scales (Shen et al. 2022). Compared with ground-observed phenological data, satellite phenological data typically have coarser spatial resolutions (ca. 0.25–8 km) and larger uncertainties brought by data sources, processing and inversion methods (Chen et al. 2015b; Wang et al. 2014). Model evaluations using ground-observed phenological data could reflect the important effects of photoperiod on herbaceous plant photomorphogenesis and foliar development at local scales. However, at regional scales, photoperiod and maximum temperature have much smaller spatiotemporal variability than precipitation and minimum temperature (Supplementary Fig. S5). Moreover, there is universal warming taking place on the QTP, so precipitation is more closely related to the changes in vegetation greenness and cover, especially in meadow and steppe regions (Ding et al. 2007; Wang et al. 2020b; Zhu et al. 2023). As a result, satellite-based model evaluations showed more prominent and extensive effects of precipitation and daily minimum temperature, while underestimated the important effects of daily maximum temperature and photoperiod.
Previous studies on the QTP have reported greater temporal variations in satellite-derived SOS than ground-observed GUD and smaller variations in satellite-derived EOS than ground-observed LSD (Shen et al. 2015b, 2022; Sun et al. 2020). The overall inter-annual change trends in satellite-derived EOS of alpine grasslands on the QTP were commonly small and insignificant (Che et al. 2014; Shen et al. 2022; Zhang et al. 2018; Zu et al. 2018). Comparison of our ground-observed and satellite-derived phenological data at several alpine grassland sites can also support this point of view (Supplementary Fig. S6). Sun et al. (2023) have reported that autumn LSD of herbaceous plants on the QTP were driven by multiple factors and the combined contribution of preseason precipitation and air temperatures to the inter-annual variability in LSD was only 38.0%. Hence, satellite-based evaluations showed that climatic factors could effectively regulate leaf onset dates, but had no significant regulating effect on leaf offset dates. Although ground-based evaluation revealed better effects of the four factors on predictions of leaf offset dates, their impacts were not as prominent as those on leaf onset dates.
Importance of precipitation in modeling foliar phenology
We used PA instead of VPD to indicate the water limitation for vegetation growth and canopy development in the AGSI model. In the original GSI, VPD is used as a surrogate for water stress because it is a continuous variable that affects soil water balance and plant stomatal closure. Sun et al. (2018) have compared the correlation between PA and VPD with soil water content as well as the GSI model performance using PA and VPD at several alpine meadow sites. As a result, clear advantages of PA over VPD were shown. To further demonstrate the effectiveness of using PA instead of VPD in modeling foliar phenology of QTP alpine grasslands, we constructed a new AGSI formulation consisting of Tmin, Tmax, VPD, and Photo (i.e. using VPD to replace PA in the AGSI). This new formulation was also calibrated and validated using the methods described in Section Model calibration and validation based on ground observations at the same four sites. The average RMSEs of predicted leaf onset and offset dates from this new formulation were 6.7 days and 7.2 days, respectively, which were higher than 6.3 days and 6.6 days from the AGSI model including PA instead of VPD.
Ground-based evaluations showed that after removing PA from the AGSI model, RMSEs of predicted leaf onset and offset dates increased by 7.7 days (145.3%) and 2.4 days (48.0%), respectively. Satellite-based evaluations found that the RMSEs of predicted leaf onset dates using MODIS and AVHRR phenological data increased by 45.3 days and 30.9 days, respectively (Fig. 7). These model-based results are consistent with previous studies based on statistical analyses that have revealed an important role of precipitation in regulating spring and autumn phenology on the QTP (An et al. 2020; Shen et al. 2015a; Sun et al. 2023). Field observations and experiments have reported that increased precipitation led to earlier green-up and later senescence in alpine grasslands regardless of warming, and the influence of precipitation was more prominent on green-up dates (Li et al. 2016; Ma et al. 2023), which is also consistent with our evaluation results. Some studies also revealed the more extensive effects of precipitation than temperatures in regulating foliar phenology of alpine grasslands on the QTP (An et al. 2020; Sun et al. 2023). This is in accord with our investigations of predominant climatic controls over the entire alpine grasslands on the plateau.
From the literature, we find that the regulatory mechanism of precipitation mainly includes two aspects. In one aspect, precipitation can substantially regulate the water condition for vegetation growth and development. Precipitation is closely correlated with soil moisture and VPD (Sun et al. 2018). Decreased precipitation may cause water deficit and stress, which could inhibit plant germination and photosynthesis activities, accelerate chlorophyll and protein degradation, and lead to premature senescence of leaves (Ren et al. 2019). Therefore, at large spatial scales, precipitation controls the inter-annual and inter-site changes in vegetation productivity and canopy development, particularly for alpine grasslands that have shallower roots than other grasslands (Chen et al. 2015b; Zhai et al. 2022). The other aspect is that precipitation generally plays a critical role in the temperature dependency of vegetation phenology in arid and semiarid regions (Shen et al. 2015a). On the QTP, it is likely to have a significant positive influence on the sensitivity of foliar phenology to air temperature changes (Shen et al. 2022; Sun et al. 2023).
Different effects of daily maximum and minimum temperatures
Both daily maximum and minimum temperatures were used in the AGSI model considering the different effects of maximum and minimum temperatures on foliar phenology. Previous models commonly used daily mean or minimum air temperature, and all existing GSI models have used daily minimum temperature. Our evaluations based on the AGSI model showed that the minimum temperature was indeed a critical factor determining foliar phenology of alpine grasslands on the QTP. The maximum temperature was also an important limiting factor because excluding its effects from the AGSI model, RMSEs of predicted leaf onset and offset dates increased by 2.7 days (50.9%) and 2.8 days (56.0%), respectively. We also compared performance of the AGSI model considering both effects of daily minimum and maximum temperatures with those revised formulations incorporating only the minimum or maximum temperature. The AGSI model showed the lowest RMSEs regardless of phenological data used, which demonstrates that incorporating effects of both minimum and maximum temperatures into a model improved prediction accuracy of foliar phenology.
Our evaluations based on both ground and satellite data showed an overall weaker limitation of air temperature than precipitation on foliar development of QTP alpine grasslands. This is contrary to a lot of previous studies that have demonstrated the dominant temperature control on vegetation phenology (Cong et al. 2017; Ding et al. 2016; Qin et al. 2022; Zhang et al. 2018; Zhu et al. 2017), but consistent with quite a number of studies that have reported a more important and extensive effect of preseason precipitation (An et al. 2020; Li et al. 2016; Ma et al. 2023; Shen et al. 2015a; Sun et al. 2023). Our model-based research also showed that compared with the maximum temperature, the minimum temperature had a more prominent influence on accurate predictions of foliar phenology, which corresponds with most studies on the QTP (An et al. 2020, 2022; Lang et al. 2019; Yang et al. 2017). However, it does not mean that daily maximum temperature is not important in the phenological prediction. Model-based evaluations have demonstrated the effectiveness of incorporating maximum temperature and revealed different effects of maximum and minimum temperatures on predicting leaf onset and offset dates.
Based on the AGSI model, we found that daily minimum temperature exerted a greater influence on leaf onset dates, while the maximum temperature showed a more important effect on leaf offset dates. Previous studies have revealed that plants need to meet the chilling and heat requirements to germinate or green-up (Chuine 2000; Wang et al. 2020a; Zhang et al. 2022). Daily minimum temperature is closely related to the chilling accumulation (Meng et al. 2019), and underlying eco-physiological processes before germination, such as transformation from starch to sugars, are more affected by the minimum temperature (Zhang et al. 2022). Higher minimum temperature can not only reduce the occurrence of freezing damage in early spring, but also mitigate the low-temperature constraints on heat accumulation before leaf onset (Shen et al. 2016). Hence, daily minimum temperature has substantial impacts on spring foliar phenology. Maybe both spring and autumn foliar phenology are more sensitive to the minimum temperature change, but there are too many constraints on autumn leaf senescence (Sun et al. 2023), which leads to the more distinct effects of minimum temperature on leaf onset dates. The constraints of daily maximum temperature were prominent in summer, and thus had a considerable influence on autumn leaf senescence (Supplementary Fig. S3). On the one hand, daily maximum (daytime) temperature controls the rates of photosynthesis and transpiration, so very low daily maximum temperature restricts vegetation production and canopy development. This effect is more popular in alpine meadows where precipitation is more and soil is wetter. On the other hand, very high maximum temperature can lead to low topsoil moisture and decreased water availability. It causes arid and semiarid regions (e.g. alpine steppes) are prone to water shortage and drought, which is also an important reason for leaf senescence and fall (Wu et al. 2022; Zhu et al. 2023). Calculation of iTmax in the AGSI considers both effects of insufficient and excessive maximum temperature. We once tried a simple linear function, like that of iTmin, to represent the limiting effect of maximum temperature, and clearly lower modeling accuracy was found (Supplementary Table S7), indicating that a relatively complex representation of the effects of maximum temperature was not unnecessary.
Limitations of the model
Although parameters in the AGSI model have specific meanings, we could only acquire those parameter values through mathematical optimization since the true values are seldom measured. The optimized values are expected to be close to the true values, yet the insufficiency of observational phenological data may lead to biases in model calibration and validation processes. The optimized parameter values are also related to phenological data used. We compared the parameter values acquired from ground-observed and satellite-derived phenological data, and found distinct differences in some parameters (Supplementary Table S5). This is the main reason causing different patterns of predominant climatic controls between ground-based and satellite-based evaluations. In addition, same model parameters were used to predict both leaf onset and offset dates under the GSI modeling framework. We once conducted an experiment to explore the difference in model performance between using one and using two AGSI thresholds (Supplementary Fig. S7). No significant difference was found, so we used a single AGSI threshold to decrease the number of model parameters.
The AGSI model only considered four climatic factors and their linear influences on foliar phenology, which also brings uncertainties to this study. Although the four factors we included in the AGSI model are all important controls of foliar phenology that have been demonstrated by many previous studies, they could not totally explain variations in leaf onset and offset dates, and prediction errors are over six days. This may be because foliar phenology of alpine grasslands is also influenced by other factors, such as preceding phenophases, soil water content, grassland degradation, community succession, nitrogen addition, and grazing (Peng et al. 2021; Shen et al. 2022). It should be noted that the predominant climatic controls we investigated are just the most limiting factors among precipitation, maximum air temperature, and minimum air temperature, not considering other factors. The limiting effect of each factor was depicted by a linear or piecewise linear function. In fact, the climatic constraints probably follow non-linear curves, and models may perform better by addressing the non-linear relations. Therefore, further efforts are needed to accurately address the effects of primary climatic factors.
There are also limitations in the model evaluation and application that need to be clarified. The AGSI model is not a process-based phenological model, so we only compared the model performance with those models of the GSI family that had been applied to grasslands worldwide. Although we did not include process-based models in our model performance evaluation, the accuracy of the AGSI model and process-based phenological models is overall comparable when assessed using ground observations (Chen et al. 2015b; Lang et al. 2019). However, process-based models developed on the QTP are generally used at the species level, which is different from the application scales of the AGSI model. Like other GSI models, the improved AGSI model can be widely applied at site-to-global scales, not at the species level, because AGSI represents an integrated climatic constraint on vegetation canopy development (Jolly et al. 2005; Sun et al. 2018). Strictly speaking, the AGSI model should not be compared with those species data-driven process-based models. Because of this reason, we used satellite phenological data to evaluate model performance at more sites in this study. It is better to use satellite phenological data obtained specifically on the QTP. However, phenological data inversion is beyond the scope of this study focusing on phenological modeling. So, we used two easily accessible and widely used satellite phenological products directly, though there may be some uncertainties in these satellite phenological data because their inversion processes are not specifically designed for the QTP (Chen et al. 2015b; Sun et al. 2020).
Conclusions
This study proposed an improved GSI model, which was referred to as the AGSI model, for predicting foliar phenology of alpine grasslands. The model was driven by daily minimum air temperature (Tmin), daily maximum air temperature (Tmax), precipitation (PA), and photoperiod (Photo). Evaluations of model performance using ground and satellite phenological data showed that the AGSI model had higher accuracy in predicting leaf onset and offset dates. Based on the AGSI model, we further assessed the impacts of Tmin, Tmax, PA, and Photo on modeling accuracy, and identified the predominant climatic controls over the entire alpine grasslands on the QTP. All four factors were found to have important effects, though there were distinctions in the results evaluated using ground and satellite phenological data. Overall, PA and Tmin showed more critical and extensive constraints on foliar development, and the limiting effect of Tmax was also considerable, especially during July–November. This study not only provides an effective model for predicting foliar phenology of alpine grasslands, but also offers a simple tool for evaluating the climatic effects on vegetation canopy development.
Supplementary Material
Supplementary material is available at Journal of Plant Ecology online.
Table S1: Geographical coordinates, elevations, mean thermal–moisture conditions, and grassland types at the 32 alpine grassland sites on the Qinghai–Tibetan Plateau (QTP).
Table S2: Grassland types, observed plant species, and time spans of observation data at the four ground sites on the QTP.
Table S3: Parameters in the Alpine Growing Season Index (AGSI) model.
Table S4: Different formulations that were constructed for assessment of the impacts of individual factors in the AGSI on modeling accuracy.
Table S5: Optimum parameter values used in the investigation of predominant climatic controls on foliar phenology of entire alpine grasslands on the QTP.
Table S6: The impacts of individual factors in the AGSI on modeling accuracy, which were average results calculated using ground-observed phenological data at four ground sites (i.e., Qumarleb, Henan, Gade, and Haiyan) and optimum model parameter values minimized total RMSEs of predicted leaf onset and offset dates.
Table S7: Comparison of modeling accuracy of the AGSI (i.e., the new model proposed in this study) and AGSInew (i.e., a new AGSI using a simpler calculation formula for iTmax), which were evaluated using ground observation data at four sites.
Figure S1: A flow chart showing the data and methods used in this study to evaluate the model performance and robustness.
Figure S2: An experiment aimed to find the optimum precipitation accumulation period for calculating PA.
Figure S3: The multi-year average daily index values of iTmax, iTmin, iPhoto, iPA, iAGSI, and AGSI at (a) Qumarleb, (b) Gade, (c) Henan, and (d) Haiyan, showing the seasonal climatic constraints.
Figure S4: Comparison of (a, c, e, g) simulated leaf onset dates and observed green-up dates, and (b, d, f, h) simulated leaf offset dates and observed leaf senescence dates at Qumarleb, Gade, Henan, and Haiyan sites.
Figure S5: Coefficients of variation (CV) of the four factors used in the AGSI model, including daily maximum air temperature (Tmax), minimum air temperature (Tmin), precipitation (Pre), and photoperiod (Photo).
Figure S6: Comparison of ground-observed and satellite-derived foliar phenology at (a) Qumarleb, (b) Gade, (c) Henan, and (d) Haiyan sites.
Figure S7: Comparison of the model accuracy between using one and using two AGSI thresholds to identify leaf onset and offset dates.
Authors' Contributions
Qingling Sun (Conceptualization, Formal analysis, Funding acquisition, Methodology, Visualization, Writing—original draft, Writing—review & editing), Jiang Zhu (Investigation, Methodology, Visualization), Baolin Li (Conceptualization, Methodology, Writing—review & editing), Siyu Zhu (Investigation, Visualization), Jie Zhu (Data curation, Resources, Writing—review & editing), Xiuzhi Chen (Resources, Writing—review & editing), and Wenping Yuan (Funding acquisition, Supervision)
Funding
This work was jointly funded by the National Natural Science Foundation of China (42201059), the General Program of Guangdong Provincial Natural Science Foundation (2024A1515012731), and the Science and Technology Program of Guangdong (2024B1212070012).
Conflict of interest statement. The authors declare that they have no conflict of interest.