# Empirical Approach for Modelling Tree Phenology in Mixed Forests Using Remote Sensing
## 2 Materials and Methods
### Data Analysis and Modelling  
The 8-day aggregate GPP data were used to reconstruct continuous daily GPP values. First, a fraction of 1/8 was applied to obtain the 8-day average from the downloaded 8-day sums. Then, each year of data was fitted with the best (low RMSE) of (i) a Gaussian function, (ii) a double logistic function after [25], and (iii) another double logistic function after [26]. (ii) and (iii) were performed using the greenbrown package (version 2.5.0) in R (see function equation details in the respective studies), while (i) was performed using non-linear least squares regression also in R, by fitting the following standard Gaussian function:  
\[\mathrm{GPP}=\mathrm{k}*\exp\left(-\frac{(\mathrm{t}-\mu)^{2}}{2\sigma^{2}} \right). \tag{1}\]  
where k is the peak GPP for the year, t is the day of the year, \(\mu\) is the day of the year corresponding to the peak of GPP, and \(\sigma\) is the standard deviation.  
The method used to extract SOS and EOS from the generated daily data is based on the widely used seasonal midpoint (SM) method (illustration in Figure 3). The traditional SM method derives SOS and EOS as the day of the year corresponding to the midpoint (0.5 * amplitude), i.e., SOS on the ascending part of the curve, and EOS on the descending part. For more stable estimates, several studies used a multi-year median of midpoints, instead of a different midpoint for each year [12; 37].  
In this study, we also used the long-term median threshold. In addition, instead of the traditional midpoint threshold (0.5), we allowed the threshold to vary from 0.05 to 0.95, and to be different for SOS and EOS. This was based on the assumption that in a mixed forest, the species of interest may start contributing significantly to the increase or decrease in GPP at different thresholds compared to a pure forest, depending on the coverage of the species of interest and also the other species present on the site. The best threshold is the one that gives the lowest RMSE of the observed versus modelled SOS or EOS. To test the predictive abilities of the approach, a leave-one-out cross-validation was performed,  
Figure 3: Example of one-year series of MODIS GPP data (_Fagus sylvatica_ at Bohinjska for year 2015). Grey dots are the 8-day average for the MODIS GPP 500 m pixel, the red line is the fitted line, a is the amplitude, SOS and EOS represent the start of season and end of season, respectively, determined using the traditional SM method.  
i.e., estimating each time the best threshold based on all years but one, and then using the holdout year as test data. In that way, all years are used exactly once as validation data when they are not used to train the model. All predictions of holdout year are then put together to calculate RMSE and graphically represent the performance of the model, hence giving a good idea of how well the approach can work in future years.  
After the cross-validation, all years were used to determine the best threshold to be applied for future years, for each combination of stations and species. To test the effect of land cover on the best thresholds obtained, fractional land share (% of the pixel) of the species of interest, other broadleaved species, coniferous species, and other land use (including agriculture, grassland, bare land, buildings, etc.) was computed. The best thresholds of 0.45, 0.5 and 0.55 (i.e., -the traditional threshold) were compared to other best thresholds (<0.45 or >0.55), using the Mann-Whitney test [50; 51].