Plant Phenology Simulation and Trigger Threshold Based on Total Climatic Production Factors – a Case Study of Stipa krylovii Phenology

1. Introduction

Plant phenology is a mutually adaptive growth and development rhythm formed by long-term adaptation of plants to seasonal changes in environmental conditions such as temperature, precipitation, and light [1], which not only reflects plant growth and development but also indicates climate change [2] and has been widely used to guide agricultural activities and disaster prevention and mitigation [3]. Meanwhile, plant phenology affects water-heat exchange and carbon cycling in ecosystems [4,5], and is an important parameter in land surface process models and plant productivity models [6]. Therefore, it is important to study plant phenology patterns and their relationship with environmental conditions.

It was found that plant spring phenology tends to be earlier and autumn phenology tends to be later in a warming context [7,8], and this phenomenon was more pronounced in the high latitudes of the Northern Hemisphere [9]. However, there is a clear spatial variability in the response of grassland plant phenology to climate change in Inner Mongolia, with a trend of earlier green-up and later wilting in the southern region, and a trend of later green-up and earlier wilting in the central and northern regions [10], and the factors causing spatial variability in grassland plant phenology are still unclear. Meteorological factors that affect plant phenology include temperature, precipitation, light, air humidity, carbon dioxide concentration, etc. Among them, temperature is considered to be the most important environmental factor affecting plant phenology [11], and plants can only grow and develop in a certain temperature environment and need a certain cumulative temperature to complete their life cycle [12]. Water deficit limits the use of light and heat conditions by plants [13] and is considered a key factor in regulating vegetation activity in arid and semi-arid regions [14]. Light is the source of energy for photosynthesis in plants, and the organic matter produced by photosynthesis is the material basis for plant growth and development [15], and photoperiod has also been shown to be an important indicator of how light affects plant phenology [16]. Meanwhile, there are significant interactions between temperature, precipitation, and light, which together affect the changes in plant phenology [17,18]. Phenology models have evolved from statistical to mechanistic models, and have been widely used for simulation of phenology and prediction of future climate change impacts in areas lacking phenology observations [19]. However, limited by the understanding of the mechanisms influencing phenological changes, especially because most of the existing studies only consider the phenological effects from single or several environmental factors on specific plant species, the response of plant phenology to the combined effects of environmental factors and its mechanisms are still unclear, making these models are not able to reflect the realism of vegetation growth and cannot effectively simulate phenological periods [20].

The research showed that plant phenology was closely related to the dynamics of photosynthesis [21]. Plant photosynthesis is the result of the interaction between environmental factors and plant biological properties, reflecting the influence of the total climate production factors. Climate is the most important factor affecting plant growth and development and is the basis for morphological establishment and physiological and biochemical changes in plants [22,23]. Climatic production potential refers to the highest biological or agricultural yield per unit area of land when other conditions (e.g., soil, nutrients, carbon dioxide, etc.) are at optimum conditions and the local climatic resources such as light, heat, and water are fully and rationally utilized [24]. Climate production potential not only reflects the influence of climate factors (temperature, moisture, light, etc.) and their combined effects on plant production, but also ensures the uniformity of the influencing factors during the whole process of plant growth and the cyclical changes in their interaction with the environment; at the same time, climate production potential also reflects the combined effects of biological factors (e.g., leaf area), environmental factors and their interactions, and can reflect the effects of extreme weather and climate events. Therefore, the use of climate production potential as a driver of plant phenology changes can avoid the shortcomings of existing models and achieve accurate simulation of phenology [25].

Chinese temperate grasslands are the third largest in the world [26], sensitive to climate change, and play an important role in the global carbon cycle [27]. The grassland of S. krylovii is one of the representative types of typical grasslands [28], which occupies an important position in livestock production [29] and has been significantly affected by warm and dry climate [30]. Using the long-term phenological and corresponding meteorological observation data from 1985 to 2018 in the grassland of S. krylovii, this study proposes the new concepts of climate production potential, cumulative climate production potential (reflecting resource accumulation), first-order derivative of cumulative climate production potential (reflecting the rate of resource change) and second-order derivative of cumulative climate production potential (reflecting sudden resource change) based on the total climate production factors influencing plant phenological changes [25]. On this basis, this study will intend to (1) verify the feasibility of simulating plant phenology and triggering thresholds based on total climatic production factors, and (2) clarify the relationship between the main phenological periods of S. krylovii and total climatic production factors and their triggering thresholds to improve the understanding of the response of phenology to the combined effect of meteorological conditions and provide a basis for the development of phenological models.

2. Materials and Methods

2.1. Study Area and Data

The study data are from the Xilinhot National Climate Observatory in Xilinhot, Inner Mongolia, China (44^◦08′ 03′ N, 116^◦19′ 43′ E, 990 m a.s.l.). It is located in the middle of the Inner Mongolian Autonomous Region, which is a typical temperate semi-arid continental climate zone. The mean annual temperature and mean annual precipitation were 2.0℃and 260 mm, respectively. Winter is cold and dry, summer is warm and humid, and solar radiation is strong. The soil type was chestnut soil. The dominant species of the ecosystem is S. krylovii, accompanied by Leymus chinensis (Trin.) Tzvel.), Allium tenuissimum Linn., Cleistogenes squarrosa (Trin.) Keng, Artemisia frigida Willd., Allium amsopodium Ledeb., Kochia prostrata (Linn.) Schrad., Artemisia scoparia Waldst. Et Kit., Heteropappus altaicus (Willd.) Novopokr. The phenology and meteorological data have been observed at the Xilinhot National Climate Observatory, Inner Mongolian, China Meteorological Administration since 1985. All meteorological data were examined and verified by the National Meteorological Information Center of the China Meteorological Administration. Plant phenology was obtained at the natural pasture observation site from the Xilinhot National Climate Observatory, Inner Mongolian, China Meteorological Administration. The fence area of the observation field is 100 m×100 m, which is divided into 4 plots of 50 m×50 m, and each plot is divided into 4 replicates for observation. Phenological periods, which were recorded by local professional observers in accordance with standard observations [31], were observed every 2 d. And 10 plants with respectable growth and complete life history for 3 consecutive years were selected for observation. A phenology is determined by the date on which 50% of the observed plants reach the phenology. In this study, the phenological dataset of dominant species S. krylovii in typical grassland include the green-up date, the heading date and the wilting date from 1985 to 2018. The phenological observation records are converted into the daily sequence from January 1st by using the method of Julian day. The meteorological dataset includes daily mean air temperature (℃), precipitation (mm), sunshine hours (h), average pressure (hpa), mean wind speed (m/s) at the height of 10m and relative humidity (%) from 1985 to 2018.

2.2. Method

2.2.1. Climate tendency rate

The climate tendency rate reflects the changing trend of meteorological elements in a region. X_i is a meteorological element value at the Xilinhot National Climate Observatory from 1985 to 2018, T_i is the corresponding chronological order, and the linear regression equation is as follows:

$$X_i=a+b{T}_i$$

(1)

where a is the regression constant, b is the regression coefficient, and b × 10 is the climate tendency rate. If the climate tendency rate is positive, it indicates that the element is showing an increasing trend; if it is negative, it indicates a decreasing trend.

2.2.2. Climate production potential

Climate production potential (LNPP) refers to the climate production per unit leaf area [25], which is obtained by dividing the climate production (YW) by the leaf area correction function f(L), thus making the climate production potential comparable among different phenological periods.

$$LNPP=\raisebox{1ex}{$YW$}\!\left/ \!\raisebox{-1ex}{$f(L)$}\right.$$

(2)

The climate production (YW) is calculated based on the step-by-step correction method recommended by the Food and Agriculture Organization of the United Nations (FAO) [25]. Firstly, the photosynthetic production (YQ) is estimated based on plant physiological mechanisms and energy conversion. Secondly, the light temperature production (YT) is evaluated by temperature corrections on YQ. Finally, the climate production (YW) is obtained by moisture corrections on YT. YW can be calculated based on daily scale meteorological data during plant growth:

$$YW=Q\times f\left(Q\right)\times f\left(t\right)\times f\left(w\right)=YQ\times f\left(t\right)\times f\left(w\right)=YT\times f\left(w\right)$$

(3)

where YW is the climate production per unit area (kg hm⁻²); Q is the total solar radiation during the plant growth period (MJ m⁻²); f(Q) is photosynthetic efficiency coefficient (g kJ⁻¹), f(t) is the temperature correction coefficient, f(w) is the moisture correction coefficient; YQ is the photosynthetic production potential per unit area (kg hm⁻²); YT is the production potential of light and temperature per unit area (kg hm⁻²).

YQ is the plant yield determined solely by solar radiation under assumed optimal conditions such as temperature, moisture, soil conditions, CO₂ concentration, and agricultural facilities, calculated as follows:

$$\begin{array}{l}YQ=f\left(Q\right)\times Q\\ =C\times S\times \epsilon \times \phi \times \left(1-\alpha \right)\times \left(1-\beta \right)\times \left(1-\rho \right)\times \left(1-\gamma \right)\\ \times \left(1-\omega \right)\times f\left(L\right)\times E\times {\displaystyle \sum Q_i\times {\left(1-\eta \right)}^{-1}\times {\left(1-\xi \right)}^{-1}}\times q^{-1}\end{array}$$

(4)

where Q is the total solar radiation during the plant growth season (∑Q_i) (MJ m⁻²), f(L) is the leaf area correction function. The other parameters listed in the equation are constants related to the biological characteristics of the plant itself (Table 1).

Daily solar radiation ($Q_i$) could be expressed as [32]:

$$Q_i=Q_0(a+b\frac{n}{N})$$

(5)

where $Q_0$ is the astronomical radiation (MJ·m^–2), n and N are the actual daily sunshine hours (h) and the maximum possible sunshine hours (h), respectively. a and b are empirical coefficients, taken as 0.29 and 0.557 [33].

$$Q_0=\left(T\times \frac{I_0}{\mathsf{\pi }{\rho }^2}\right)\times \left(\omega \sin \phi \sin \delta +cos\phi cos\delta \sin \omega \right)$$

(6)

$$N=24\times \frac{\omega }{\pi }$$

(7)

where T is the time in days (=24×60 min), I₀ is the solar constant (0.082 MJ·m⁻²·min⁻¹), ω is the time angle (rad), φ is the geographic latitude (rad), δ is the solar declination (rad), ρ is the distance between the sun and the earth, ρ² is the correction coefficient of eccentricity of earth orbit.

$$\frac{1}{{\rho }^2}=1.00011+0.342221\cos \theta +0.00128\sin \theta +0.000719\cos 2\theta +0.00077\sin 2\theta$$

(8)

$$\omega ={\cos }^{-1}[-\tan \phi \tan \delta ]$$

(9)

$$\begin{array}{l}\delta =0.006918-0.399912\cos \theta +0.070257\sin \theta -0.006758\cos 2\theta \\ +0.000907\sin 2\theta -0.002697\cos 3\theta +0.000148\sin 3\theta \end{array}$$

(10)

$$\theta =\raisebox{1ex}{$2\pi (d_n-1)$}\!\left/ \!\raisebox{-1ex}{$365$}\right.$$

(11)

where θ is the solar angle (rad), d_n is the chronological order, with January 1st as 1, January 2nd as 2, and so on.

The light temperature production (YT) is the plant yield determined solely by solar radiation and air temperature under assumed optimal conditions such as moisture, soil conditions, CO₂ concentration, and agricultural facilities, calculated as follows:

$$YT=YQ\times f\left(t\right)$$

(12)

$$f\left(t\right)=\{\begin{array}{l}0,t<t_{\min },t>t_{\max }\\ \frac{t-t_{\min }}{t_s-t_{\min }},t\le t<t_s\\ \frac{t_{\max }-t}{t_{\max }-t_s},t_s\le t<t_{\max }\end{array}$$

(13)

where $\left.f(t\right)$ is the temperature correction coefficient, 0≤f(t)≤1, t is the daily average temperature (℃), t_min, t_s, t_max are the lower limit temperature, optimal temperature, and upper limit temperature (℃) of the growth period, respectively, where 0 ℃, 20 ℃ and 35 ℃ are taken in this study [34].

The climate production (YW) is the plant yield determined solely by solar radiation, air temperature and moisture under assumed optimal conditions such as soil conditions, CO₂ concentration, and agricultural facilities, calculated as follows:

$$YW=YT\times f\left(w\right)$$

(14)

$$f\left(w\right)=\{\begin{array}{l}1,VPD\text{<}VP{D}_{\min }\\ \frac{VP{D}_{\max }-VPD}{VP{D}_{\max }-VP{D}_{\min }},VP{D}_{\min }\le VPD\le VP{D}_{\max }\\ 0,VPD\le VP{D}_{\max }\end{array}$$

(15)

$$VIP=0.6108\exp [\frac{17.27t}{t+237.3}](1-RH)$$

(16)

where VPD is the water vapor pressure deficit (Pa), t the daily average temperature (℃), RH is the average relative humidity (%), VPD_max and VPD_min are constants determined by vegetation type, VPD_max and VPD_min are taken as 3100 Pa and 650 Pa, respectively [35].

To achieve the prediction of plant phenology, it is necessary to calculate the LNPP for each day starting from January 1st and calculate the cumulative LNPP to reflect the climate resource status, i.e., the cumulative climate production potential. The cumulative climate production potential is a function of time and conforms to the Logistic curve. The first derivative and second derivative can be obtained according to the cumulative climate production potential, which respectively reflect the change speed and acceleration (Abrupt change) of climate resources. The mutation (second subderivative) of the cumulative climate production potential can be used to indicate the beginning and the ending of plant growth, while the cumulative climate production potential or its change speed (first subderivative) can be used to indicate the other phenological dates between the beginning and the ending of plant growth [25].

2.2.3. Logistic curve

The accumulation of daily climate production potential can reflect the degree to which plants can utilize climate resources. Thus, the relationship between plant phenology and climate production potential can be analyzed. To quantitatively describe the accumulation of climate production potential over time, the Logistic curve is used to simulate it. The trend equation of accumulative climate production potential based on the Logistic curves is calculated as follows:

$$\frac{dY}{dt}=rN(\frac{K-Y}{Y})$$

(17)

The integral equation could be expressed as:

$$Y=\frac{K}{1+a{e}^{-rt}}$$

(18)

where dY/dt is the instantaneous growth rate, $K$ is the theoretical upper limit of climate production potential, $t$is time (d), Y is the cumulative value of climate production potential from January 1st, a is the coefficient, r is the intrinsic growth rate, and is the maximum instantaneous growth rate under specific conditions.

2.3. Data processing

Microsoft Excel 2019 is used to collate and summarize the data, Matlab 2018a software is used to fit the Logistic curve and calculate the change rate (the first derivative of the cumulative climate production potential) and mutation (the second derivative of the cumulative climate production potential) of the cumulative climate production potential, and SPSS 21 is used to analyze the change trend of the main phenological dates of S. krylovii and their relationships with the cumulative climate production potential, its change rate and mutation. The Origin 2018a is used to plot the figures.

3. Results

3.1. Climate change trends

From 1985 to 2018, the annual average temperature of the Xilinhot National Climate Observatory in Inner Mongolia showed an upward trend, while annual precipitation and annual sunshine hours showed a downward trend (Figure 1). From 1985 to 2018, the average annual temperature in the region was between 1.3 and 4.8 ℃, with an average of 3.3 ℃, showing a significant upward trend (0.40 ℃/10a, P<0.01). The average temperature from December to January of the following year showed an insignificant downward trend, while the average temperature from February to November showed an upward trend. Among them, the changes in average temperature in March, July, August, and September reached a significant level, with values of 1.17 ℃/10a (P<0.01), 0.89 ℃/10a (P<0.05), 0.53 ℃/10a (P<0.05), and 0.55 ℃/10a (P<0.05), respectively. The interannual variability of annual precipitation was large, ranging from 121.1 to 511.7 mm, with an average annual precipitation of 272.4 mm, showing an insignificant decreasing trend (–8.24 mm/10a, P>0.05), mainly occurring in July and August (–9.73 mm/10a, P>0.05; –10.76 mm/10a, P>0.05). The annual sunshine hours range from 2720.8 to 3215.8 hours, with an average of 2991.9 hours, showing an insignificant downward trend (–32.10 h/10a, P>0.05). Except for January, February, July, and August, the average monthly sunshine hours in the Xilinhot area showed a decreasing trend, with significant changes in May and June (–9.06 h/10a, P<0.05; –10.76 h/10a, P<0.05).

3.2. Phenological change trends of S. krylovii plant

From 1985 to 2018, the green-up of S. krylovii plant in the Xilinhot area of Inner Mongolia mainly happened in mid to late April. The earliest green-up date took place in 1993, which happened on April 2, and the latest green-up date took place was in 2009, which happened on May 7. The green-up date showed a significant delay trend, with an average delay of 5.4 d/10a (P<0.05) (Figure 2). The interannual differences in the heading dates are significant, occurring from early July to middle September, with a difference of nearly three months between the earliest and latest heading dates. The heading date showed an insignificant trend of advancing, with an average of 5.4 d/10a (P>0.05). After September, S. krylovii plant gradually entered the wilting period, mainly concentrated in late September and early October, and the wilting date at the latest in middle October (1987, 1998, 2000). The wilting dates showed a weak trend of advance, with an average advance of 1.2 d/10a (P>0.05). Overall, the significant delay in the green-up date and the advance of the wilting date. As a result, the length of the growth season of S. krylovii plant showed in a significant reduction with an average reduction of 6.3 d/10a (P<0.05).

3.3. Relationship between the main phenology and climate production potential

The cumulative climate production potential of S. krylovii plant showed a Logistic curve over time (Table 2), and the determination coefficients of the cumulative climate production potential of S. krylovii plant from 1985 to 2018 were all greater than 0.996.

3.3.1. Relationship between the green-up date and climate production potential

The relationship between the green-up date of S. krylovii plant and the cumulative climate production potential indicates (Figure 3) that there is a significant positive correlation between the green-up date and the cumulative climate production potential. The correlation coefficients (R²) are 0.667, 0.822, and 0.862 for the cumulative climate production potential, its first derivative and its second derivative, respectively. Among them, the relationship between the green-up date and the second derivative of cumulative climate production potential was the best, indicating that the drastic changes in meteorological conditions were important reasons for triggering the green-up phenology of S. krylovii plant. Based on this relationship, it can be concluded that the triggering threshold of the green-up period depended on both resource change parameters and biological characteristic parameters, namely the slope (0.085) and intercept (–5.363) of the second derivative of cumulative climate production potential with diurnal variation.

3.3.2. Relationship between the heading date and climate production potential

There is a significant positive correlation between the heading date and the cumulative climate production potential of S. krylovii plant (Figure 4), but there is a significant negative correlation with the first derivative and the second derivative of cumulative climate production potential. The correlation between the heading date and the cumulative climate production potential of its first derivative and its second derivative (0.840) was significantly better than the correlation with its first derivative and second derivative (0.754 and 0.607), reflecting that the heading period of S. krylovii plant mainly depended on the cumulative degree of climate resources. Based on the relationship between the heading date and the cumulative climate production potential, it can be concluded that the triggering threshold of the heading date depended on both resource change parameters and biological characteristic parameters, namely the slope (394.632) and intercept (–38026.268) of the cumulative climate production potential with diurnal changes.

3.3.3. Relationship between the wilting date and climate production potential

There is a significant positive correlation between the wilting date of S. krylovii plant and the cumulative climate production potential and its second derivative, but a significant negative correlation with the first derivative of cumulative climate production potential (Figure 5). The correlation coefficients between the wilting date and the cumulative climate production potential, its first derivative and second derivative were 0.315, 0.917, and 0.915, respectively. The correlations between the wilting date and the first and second derivative of cumulative climate production potential were good, however, but the wilting date of S. krylovii plant showed a significant negative correlation with the first derivative of cumulative climate production potential, indicating the decrease of utilizing climate resources rather than the stop. Considering the variability between the green-up phenology and the wilting phenology [25], the sudden change in cumulative climate production potential can better reflect the changes in plant withering and yellowing periods. Based on the relationship between the wilting date of S. krylovii plant and the second derivative of the cumulative climate production potential, it can be concluded that the triggering threshold of the wilting date also depended on both resource change parameters and biological characteristic parameters, namely the slope (0.086) and intercept (–27.620) of the second derivative of cumulative climate production potential with diurnal variation.

4. Discussion and Conclusion

Meteorological factors are important factors that affect plant growth and development, as well as the basis for plant morphogenesis, physiological and biochemical changes. The meteorological factors affecting the green-up period of grassland plants include temperature, sunshine duration [36], accumulated temperature [19,37] and moisture [38,39]. When the air temperature is below the threshold, the autumn phenology occurs [40]; when the photoperiod is shortened to the threshold that limits plant growth and development, it will induce plant leaf senescence and enter a dormant state. The in situ simulation experiment of S. krylovii grassland also showed that it is environmental factors rather than plant productivity that drives leaf senescence [42]. In this study, the abrupt relationship between the green-up date and wilting date of S. krylovii plant and the second derivative of cumulative climate production potential was the best, indicating that different phenological periods of grassland plants were comprehensively affected by meteorological conditions, and sudden changes in environmental factors trigger the green-up and wilting of plants. The heading date of S. krylovii plants had the best relationship with the cumulative climate production potential, indirectly proving the impact of rising summer temperature as a driver of climate production potential on the phenological periods between the beginning and ending phenology, reflecting the cumulative climate resource effects.

In this study, the long-term phenology and corresponding meteorological observation data of S. krylovii grassland from 1985 to 2018 were used to verify the indication of climate production potential based on total climate production factors on plant phenology, and the relationship between the main phenological stages of S. krylovii plant and climate production potential was clarified. The climate production potential not only reflects the comprehensive effect of environmental factors in plant production and their consistency throughout the entire process of plant growth, but also reflects the effects of biological factors, environmental factors and their interactions, as well as extreme weather and climate events, which will be helpful for achieving accurate simulation of plant phenology. The main conclusions are as following:

(1) From 1985 to 2018, the green-up date of S. krylovii plant in Xilinhot, Inner Mongolia was delayed by an average of 5.4 d/10a, the heading date was advanced by an average of 5.4 d/10a, the wilting date was advanced by an average of 1.2 d/10a, and the length of the growing season was shortened by an average of 6.3 d/10a. The phenological change trend of S. krylovii is consistent with the existing research in this area [36,43].

(2) The climate production potential based on total climate production factors is a good indicator for plant phenology, and the sudden change in cumulative climate production potential reflects the drastic changes in climate conditions, which can effectively indicate the green-up and wilting periods of S. krylovii plant; while the cumulative climate production potential reflects the rate of climate resource change, and can reflects resource utilization to a certain degree, which can effectively indicate the heading period of S. krylovii plant.

(3) The relationship between the phenological periods and the cumulative climate production potential of S. krylovii plant further indicates that different phenological periods of grassland plants are comprehensively influenced by meteorological conditions. To comprehensively understand the response of plant phenology to environmental changes, it is necessary to understand the influencing mechanism of total climate production factor on plant phenology.