The timings of key physiological events in the Earth’s terrestrial biosphere are crucial to controlling the relationship of the carbon, energy and water cycles between the land surface and the atmosphere. By understanding the timings and spatial distributions of different phenological events such as the start of the growing season (SOS), length of the growing season (LOS) and the end of the growing season (EOS), the impacts of climate change can be better understood and assessed. Previously, studies have shown temperature to be a driving factor behind changes in northern latitude phenology, whilst precipitation drives tropical regions. As global temperatures continue to increase, this is expected to lead to a growth in vegetation productivity in the Northern Hemisphere due to an earlier SOS and an increased LOS.
Time series of vegetation indices derived from satellite imagery have been used to monitor changes in the phenological development of terrestrial vegetation over the past four decades at a number of spatial scales. Indices such as the normalized difference vegetation index (NDVI), have a strong correlation with green vegetation biomass, but lack physical meaning, suffer from sensitivity to non-canopy factors, and saturate at higher values.
The Global Climate Observing System (GCOS) defined the Fraction of Absorbed Photosynthetically Active Radiation (FAPAR) as one of the fifty Essential Climate Variables (ECVs), as it is a key biological quantity that contributes to the characterization of Earth’s climate. Using an ECV to derive and monitor the phenology of terrestrial vegetation has distinct advantages over vegetation indices, since FAPAR is a bio-physical property of vegetation describing the presence and state of the photosynthetic activity. As such, we processed FAPAR products from Envisat MERIS and Copernicus Sentinel-3 OLCI covering the period 2002 – 2020 to obtain phenological metrics.
The main objectives are to:
(1) Assess the feasibility of using satellite derived FAPAR time series for deriving phenology over the Northern Hemisphere;
(2) Discuss the spatiotemporal patterns and changes in three phenological metrics: SOS, LOS and EOS,
(3) Validate the performance of the satellite products using in-situ FAPAR measurements obtained from the national ecological observatory network.

The European Environment Agency (EEA) recently presented the high-resolution vegetation phenology and productivity (HR-VPP) product group as part of the Copernicus Land Monitoring Service. These datasets serve the European Union with 10-day seasonal data for monitoring of terrestrial ecosystem dynamics, and seasonal parameters describing key phenological stages of vegetation, such as onset, end, length, amplitude, productivity, and several others. The seasonal data are derived from from 10-m Sentinel-2 observations and the novel plant phenology index (PPI), which responds to the photosynthetic production of leaf biomass. To tackle cloudiness, an improved cloud screening algorithm was implemented. Furthermore, A new time-series processing sequence was developed based on the TIMESAT software system, to generate gap-filled and smoothed seasonal trajectories, phenological parameters, and quality assurrance layers. HR-VPP also includes near real-time observations of fAPAR, LAI and NDVI.

HR-VPP presents spatial and temporal data of vegetation productivity, providing opportunities for monitoring crop yield variations, forest growth, and ecosystem productivity and health. It has also been shown to indicate vegetation responses to drought and other climate disturbances. To ensure accuracy of the final products our team tested a range of methodological approaches and datasets, and validated the data against a broad range of reference data, including PhenoCam imagery, GPP estimated from ICOS eddy-covariance data, ground-observed phenological observations of the PEP725 network, and other ancillary data. An independent validation against ground references and other satellite-based phenology products was carried out following CEOS/LPV best practices.

The HR-VPP product group provides free and user-ready continental data at the highest spatial and temporal resolution presently available for detailed assessment of local vegetation conditions across the European continent. The production currently covers the years 2017 – 2021. Outcomes are relevant to EU and international policies such as the European Green Deal Biodiversity Strategy 2030, the EU Soil Strategy, the EU Climate Adaptation Strategy and the LULUCF regulation. In an era of rapid global change it will provide an invaluable data source for innovations and science to understand and mitigate future climate and human induced vegetation changes.

High-quality estimation of land surface phenology (LSP) is becoming crucial for understanding the impacts of climate change on the ecosystem functioning. Reliable extraction of photosynthetic phenology using remotely sensed vegetation indices is important for quantifying the seasonal cycle of atmospheric CO₂ concentrations and modelling the responses of vegetation to climatic warming. However, there are still uncertainties about time series pre-processing, the accuracy and robustness of LSP methods, the impact of spatial resolution, and the interpretability of LSP metrics. We present a line of research that demonstrates the relevance of reducing the uncertainties when studying vegetation phenology responses to climatic warming.

We firstly analysed four state-of-the-art phenology methods and different time series smoothing techniques. These methods were applied to the Copernicus Global Land Service (CGLS) SPOT-VEGETATION and PROBA-V leaf area index (LAI) 1 km V2.0 time series. The resulting LSP metrics were validated with near-surface PhenoCam, eddy covariance FLUXNET data, and in situ observations of phenophases (PEP725 and USA-NPN). We found that the threshold-based method performed the best in terms of root mean square error for the start of season (30% threshold of the annual amplitude) and the end of season (40% threshold).

Based on these findings, we developed a novel threshold-based technique named maximum separation (MS), which we implemented on a cloud-based platform (Google Earth Engine). The main advantage of the MS method is that it can be directly applied to daily non-smoothed time series without any additional pre-processing steps, reducing the uncertainty associated with excessive smoothing. The implementation of the proposed method in GEE allowed for the quick processing of phenological maps produced from MODIS and Sentinel-2 at the continental scale. Such maps highlighted the necessity of decametric resolution satellite data for proper monitoring of vegetation phenology at the canopy scale in heterogeneous landscapes.

We also examined the relevance of selecting the proper spectral indices for accurate tracking forest photosynthetic phenology. We found that structural vegetation indices (NDVI and EVI) were suited for estimating the start of the photosynthetically active season in deciduous broadleaf forests, but a physiological vegetation index (chlorophyll/carotenoid index) was better for predicting the end of the photosynthetically active season in deciduous broadleaf forests and both the start and end of the season in evergreen needleleaf forests. The divergent performances were rooted in the combined control of structural and physiological regulations of carbon uptake by plants. We propose revisiting the dynamics of photosynthetic phenology using physiological vegetation indices, which has substantial implications for global plant phenology and carbon uptake investigations.

Finally, we proposed a machine learning framework for modelling phenology metrics as a function of temperature, water availability, and radiation at the global scale. Machine learning captures non-linear and non-parametric relationships between different factors that control leaf unfolding, outperforming standard processed-based and statistical models. Such framework enables the inclusion of multiple input climate variables, defined in experimental and field studies, and ranking the relevance of factors that control leaf unfolding. These results also highlight the importance of variable and LSP method selection for effectively monitoring the responses of vegetation to climatic factors.

The timing of crop development is an important information for crop modelling and crop yield forecasting. Indeed, most crops experience significant lags in space (e.g. from one country to another) and in time (from one year to another) as a function of varying agrometeorological conditions and agricultural practices. However, only a few countries monitor crop phenological stages at large scale with regular and harmonized reporting. In France, the FranceAgriMer institute publishes the weekly Céré'Obs bulletin throughout the whole agricultural season: a report of the crop development stages from sowing to harvest based on in situ observations by experts. Reported data includes estimates of the percentage of parcels for each French NUTS-3 region that have reached a given development stage. In France, the bulletin covers five crops, i.e. winter soft wheat, winter barley, durum wheat, spring barley, and grain maize. Depending on the crop, six to seven crop management and physiological stages are considered (e.g. for winter cereals: sowing, emergence, beginning of tillering, beginning of stem elongation, node 2, inflorescence, harvest). The European Union Copernicus program and the launch of the Sentinel constellations are offering a unique source of data with the capacity to perform large-scale crop monitoring at parcel level. Many authors have shown the potential of Sentinel-2 data to construct time series and to analyse their variations in space and in time. Most studies have not specifically focused their analysis on cropland. They have developed generic phenometrics, describing basic phenological stages, such as the start and the end of the season, which are rarely connected to physiological crop development stages. In this study, we propose to design and assess an automatic crop stages identification model based on Sentinel-2 time series that is calibrated with the Céré'Obs expert observations. The use of Sentinel-2 is constrained by a mandatory a-priori knowledge of the crop type in order to extract crop-specific Sentinel-2 time series. This is possible in France, using the Geospatial Aid applications (GSAA) data, publically available for past seasons.
The proposed model is based on phenology-driving variables (PDV) on which thresholds are defined to identify crop stages at parcel level. We assessed four PDV as input of the model:
a) Growing degree days (GDD) derived from weather data, not involving the Sentinel-2 data
b) Leaf area index (LAI) derived from Sentinel-2 time series
c) Green vegetation ratio, computed from the ratio of LAI as local (parcel) minimum and maximum of the LAI time series
d) Percentage of CSDM (calibrated on LAI time series) progress
Among the four PDVs, only (a) has the benefit to not depend on any a-priori knowledge, therefore can be run easily in near real time. The same applies to PVS (b), except that crop types need to be identified to extract crop specific LAI along the season. PVSs (c) and (d) require modelling of nearly the full LAI time series. The thresholds on PDV are calibrated by optimizing the matching between the distributions of the estimated stages for all parcels in a given NUTS-3 region and the distributions as provided by Céré'Obs. Kolmogorov–Smirnov test is used for the distribution matching.
For this study, we focused on soft winter wheat and grain maize over all NUTS-3 regions of France covered by the Céré'Obs report. We analysed six cropping seasons from 2016 to 2021. The Sentinel-2-based smooth time series of LAI were derived from the in-house processing of Sentinel-2 Level-2 data, involving a biophysical variable processor, aggregation at parcel level (using LPIS/GSAA data), time series outlier detection, and the Whittaker time series smoother.

Seasonal crop phenology allows understanding a crop’s metabolic cycle, its reactions to environmental influences, and its buildup of biomass, among others. Crop phenology is hence a valuable input for numerous agricultural monitoring tasks, including the assessment of climate impacts, management practices, and yield estimation. As area-wide and field-based collection of related data is not feasible, it is of great interest for the EO research community to leverage the potential of remote sensing data for deriving phenological information on cultivated crops.
Studies so far usually had reference data available for limited regions. This is sufficient for the validation of phenology estimations on fine scales. However, this leads to an insufficient error assessment when area-wide estimations are made. Thus, approaches developed for broad-scale applications often lack precise validation by using reference data or prediction units beyond the plot level (e.g, Meroni et al., 2021). Methods for estimating phenology for large areas while providing reliable plot level information are therefore needed.
We address this need by exploiting a comprehensive data set on phenological observations across Germany provided by the German Weather Service (DWD) that we match to winter wheat plots identified from the German-wide 10 m crop type map by Blickensdörfer et al. (2021). This allows to relate time series from Sentinel-1, 2, and Landsat 8 data to field observations for 700 to 800 locations for every year since 2017. We complement our remote sensing data by meteorological data provided by the DWD.
From the remote sensing imagery, we summarize the pixel values for each plot and calculate time series of the Normalized Difference Vegetation Index (NDVI), γ0 backscatter coefficient (VH & VV), and backscatter cross-ratio. We derive daily precipitation sums from DWD precipitation radar (RADOLAN) and temperature and global radiation measurements by spatially interpolating observations from DWD weather stations. To suppress noise in the remote sensing time series, we use locally weighted scatterplot smoothing (loess) and ensure a common temporal interval of the different data sources by linear interpolation. We feed the resulting time series into a one-dimensional convolutional neural network (CNN), that is trained to deliver the entry dates of four different phenological stages for winter wheat, namely stem elongation, heading, ripening, and harvest.
Preliminary results show that the chosen method is suitable for deriving phenological metrics from the combined optical, radar and meteorological time series. The time-consuming tuning of the model architecture and the selection of the most suitable input time series, however, are still pending prior to an in-depth validation. Further tests will be performed to investigate the suitability of the approach in comparison to established methods to identify strengths and limitations.
The presented approach is new to the field of crop phenology analysis in two respects: first, the use of machine learning and, in particular, deep learning to derive crop phenology has only laterally been explored, and second, the chosen model architecture allows the direct combination of optical, radar and meteorological data instead of establishing separate rule-sets and combining them ex-post. After finalizing the model, area-wide predictions will allow to analyze spatial patterns and trends of phenology in the agricultural landscapes across Germany.
References
Blickensdörfer, L., Schwieder, M., Pflugmacher, D., Nendel, C., Erasmi, S., Hostert, P., 2021. National-scale crop type maps for Germany from combined time series of Sentinel-1, Sentinel-2 and Landsat 8 data (2017, 2018 and 2019). https://doi.org/10.5281/ZENODO.5153047
Meroni, M., d’Andrimont, R., Vrieling, A., Fasbender, D., Lemoine, G., Rembold, F., Seguini, L., Verhegghen, A., 2021. Comparing land surface phenology of major European crops as derived from SAR and multispectral data of Sentinel-1 and -2. Remote Sens. Environ. 253. https://doi.org/10.1016/j.rse.2020.112232