# Characterizing Growing Season Length of Subtropical Coniferous Forests with a Phenological Model
## Materials and Methods
### Ecosystem-Scale Plant Photosynthetic Phenology Model  
This study used a GPP-based phenology model [26] to analyze the dynamic changes of vegetation phenology in each site from 2009 to 2017, and to discuss the potential impact of disturbance on the ecosystem-scale phenological processes. The GPP-based phenology model was originally developed utilizing interannual daily maximum GPP (umol CO\({}_{2}\) m\({}^{-2}\) s\({}^{-1}\)) obtained from half-hourly observations to quantify the seasonal change of canopy photosynthetic capacity [26]. Due to the large number of evergreen trees distributed in the flux source area of the sites (Table 1), however, canopy photosynthesis remains active year-round. Even though we have observed significant seasonal variation in GPP derived from 30 min data [44], models of daily maximum GPP did not yield ecologically meaningful estimates of phenology parameters. Therefore, this study applied the GPP-based phenology model using daily cumulative GPP (g C m\({}^{-2}\) d\({}^{-1}\)). In the phenological model (Equations (2) and (3)), \(A(t)\) represents daily cumulative GPP and \(k(t)\) is the GPP growth rate (daily change in GPP), as defined by the \(A(t)\) function (g C m\({}^{-2}\) d\({}^{-1}\) d\({}^{-1}\)) at day \(t\) (\(t=1,\,\ldots\,\), 365), estimated for each site by year:  
\[A(t)=y_{0}+\frac{a_{1}}{\left[1+\exp\left(-\frac{t-t_{01}}{b_{1}}\right)^{c_{ 1}}\right]}-\frac{a_{2}}{\left[1+\exp\left(-\frac{t-t_{02}}{b_{2}}\right)^{c_{ 2}}\right]} \tag{2}\]  
\[k(t)=\frac{dA(t)}{dt}=\frac{a_{1}c_{1}}{b_{1}}\frac{\exp\left(-\frac{t-t_{01}} {b_{1}}\right)}{\left[1+\exp\left(-\frac{t-t_{01}}{b_{1}}\right)\right]^{1+c_{ 1}}}-\frac{a_{2}c_{2}}{b_{2}}\frac{\exp\left(-\frac{t-t_{02}}{b_{2}}\right)}{ \left[1+\exp\left(-\frac{t-t_{02}}{b_{2}}\right)\right]^{1+c_{2}}} \tag{3}\]  
Estimated empirical fitting parameters of \(A(t)\) (\(a_{1}\), \(a_{2}\), \(b_{1}\), \(b_{2}\), \(c_{1}\), \(c_{2}\), \(t_{01}\), \(t_{02}\), and \(y_{0}\)) were then used in the \(k(t)\) function to calculate the GPP growth rate [26] by site and year. Then, these equations were used to define predicted plant green-up and senescence phases during the growing season: The Recovery Line (\(A_{RL}\)) and Senescence Line (\(A_{SL}\)) were defined as Equations (4) and (5) [26], respectively:  
\[A_{RL}(t)=k_{PRR}t+A(t_{PRD})-k_{PRR}t_{PRD} \tag{4}\]  
\[A_{SL}(t)=k_{PSR}t+A(t_{PSD})-k_{PSR}t_{PSD} \tag{5}\]  
where: \(k_{PRR}\) is the maximum GPP growth rate, which is the maximum recovery rate of GPP, \(t_{PRD}\) is the time at which \(k_{PRR}\) (day) occurs, \(k_{PSR}\) is the minimum GPP growth rate (the maximum GPP senescence rate), and \(t_{PSD}\) is the time at which \(k_{PSR}\) (day) occurs, which are calculated as follows [26]:  
\[t_{PRD}\approx t_{01}+b_{1}\ln(c_{1}) \tag{6}\]  
\[t_{PSD}\approx t_{02}+b_{2}\ln(c_{2}) \tag{7}\]  
\[k_{PRR}=k(t_{PRD}) \tag{8}\]\[k_{PSR}=k(t_{PSD}) \tag{9}\]  
From the phenology model, maximum GPP is defined as:  
\[A_{p}=\max\{\left|A(t)\,,t_{start}<t<t_{end}\right|\} \tag{10}\]  
where: \(t_{start}\) = 1, and \(t_{end}\) = 365 during non-leap years and \(t_{end}\) = 366 during leap years. SOS is then defined as the value of \(t\) at which GPP predicted by \(A_{RL}(t)\) is zero. The start of peak (SOP) of the growing season is the value of \(t\) when GPP predicted by \(A_{RL}(t)\) reaches \(A_{p}\). The end of peak (EOP) is defined by the date (\(t\)) when GPP predicted by \(A_{SL}(t)\) reaches \(A_{p}\), and EOS is the value of \(t\) at which GPP predicted by \(A_{SL}(t)\) is zero. Since the study area is located in a warmer subtropical region, the days between SOS and SOP are defined as spring, days between SOP and EOP are defined as summer, days between EOP and EOS are defined as autumn, and days between EOS and the following year's SOS are defined as winter (Figure 2) [26].