Estimating and Monitoring Land Surface Phenology in Rangelands: A Review of Progress and Challenges

Abstract

Land surface phenology (LSP) has been extensively explored from global archives of satellite observations to track and monitor the seasonality of rangeland ecosystems in response to climate change. Long term monitoring of LSP provides large potential for the evaluation of interactions and feedbacks between climate and vegetation. With a special focus on the rangeland ecosystems, the paper reviews the progress, challenges and emerging opportunities in LSP while identifying possible gaps that could be explored in future. Specifically, the paper traces the evolution of satellite sensors and interrogates their properties as well as the associated indices and algorithms in estimating and monitoring LSP in productive rangelands. Findings from the literature revealed that the spectral characteristics of the early satellite sensors such as Landsat, AVHRR and MODIS played a critical role in the development of spectral vegetation indices that have been widely used in LSP applications. The normalized difference vegetation index (NDVI) pioneered LSP investigations, and most other spectral vegetation indices were primarily developed to address the weaknesses and shortcomings of the NDVI. New indices continue to be developed based on recent sensors such as Sentinel-2 that are characterized by unique spectral signatures and fine spatial resolutions, and their successful usage is catalyzed with the development of cutting-edge algorithms for modeling the LSP profiles. In this regard, the paper has documented several LSP algorithms that are designed to provide data smoothing, gap filling and LSP metrics retrieval methods in a single environment. In the future, the development of machine learning algorithms that can effectively model and characterize the phenological cycles of vegetation would help to unlock the value of LSP information in the rangeland monitoring and management process. Precisely, deep learning presents an opportunity to further develop robust software packages such as the decomposition and analysis of time series (DATimeS) with the abundance of data processing tools and techniques that can be used to better characterize the phenological cycles of vegetation in rangeland ecosystems.

1. Introduction

Rangelands are defined as landscapes that are mainly comprised of grasslands, woodlands, shrublands and wetlands [1]. Rangelands provide forage, water and cover for grazing livestock and wildlife whose socio-economic value is substantial in many countries [2,3]. Globally, the state of rangelands has been threatened by many challenges such as biodiversity loss [4], soil erosion [5], frequent veld fires which destroy habitats [6] and most severely, the encroachment of alien invasive plants [7,8,9]. Alien invasive species typically outcompete indigenous vegetation, often replacing palatable grasses with plants that are poisonous to livestock [10]. A better understanding of the phenological dynamics of the various vegetation types in the rangelands will enable us to assess how climate variability and management practices affect various functional groups [1]. Phenology data is essential for designing and planning rangeland management systems. For instance, phenology data can be used to adjust the timing of grazing or to manage burns relative to the phenological cycles of vegetation species and for planning restoration actions, such as targeted grazing [11]. The monitoring and modeling of the changes in phenological cycles of vegetation can also help rangeland managers to make suitable and cost-effective decisions on how to adjust management strategies to optimize livestock production and other ecosystem services provided by rangelands.

The assessment of changes in phenological cycles of vegetation in rangelands can be achieved using three main types of observations: human visual surveillance using phenology networks [11,12,13], near surface measurements such as phenology cameras (PhenoCams) [1,14,15]; unmanned aerial vehicles (UAVs) or drones [16,17]; and remote sensing measurements estimated from polar orbiting and geostationary satellite sensors [18,19,20]. The monitoring of vegetation phenological cycles across regions was previously difficult to assess using ground observed phenological events due to limited spatial coverage [21]. The availability of remotely sensed data provided a long-term opportunity because it is strengthened by the ability of the data sets to present the phenological trends of vegetation at large spatial and temporal extents [22,23]. However, remote sensing does not observe vegetation only because it captures whatever is in the area covered by the field of view of the satellite sensor; thus the move to use the phrase ‘land surface phenology (LSP)’ is seen in the remote sensing literature. LSP is described as the estimation of the phenological cycles of different vegetation species using data acquired from satellite sensors [24].

The utility of several polar orbiting and geostationary satellite sensors has been explored in the estimation and monitoring of LSP in rangelands. The most common satellite data sources include Landsat [25,26], Advanced Very High Resolution Radiometer (AVHRR) [27,28], Moderate Resolution Imaging Spectroradiometer (MODIS) [19,29,30], Sentinel-2 [31], PlanetScope [1] and Himawari Imager [32]. However, the applications of each of these data sources present various challenges and opportunities in LSP research. For example, early sensors such as AVHRR and MODIS have a high temporal resolution sufficient to capture subtle changes in vegetation development, but they have low spatial resolution insufficient to capture plant specific phenological changes [33]. On the other hand, medium and high spatial resolution sensors such as SPOT and Sentinel-2 with a 3-10 day temporal resolution do not provide adequate time series data for the characterization of the phenological cycles of vegetation, especially in areas with high frequent cloud cover [34]. Recently, data fusion methods have emerged as a solution to minimize the trade-offs associated with the spatial and temporal characteristics of satellite sensors in LSP investigations [35,36,37].

In LSP investigations, scientists do not use raw spectral bands to estimate the phenological trends of vegetation; instead, they use vegetation indices (VIs) and plant biophysical variables such as the leaf area index (LAI). Vegetation indices are calculated using spectral data in the visible and near infrared (NIR) parts of the electromagnetic spectrum [38]. The spectral data in the visible and NIR sections of the electromagnetic spectrum are commonly used because they are more sensitive to plant growth and development [39]. Long term satellite data archives present an opportunity to retrospectively extract phenological characteristics of vegetation using a wide range of vegetation indices [40]. Although there is an abundance of spectral VIs that have been established for various functions, the current review focuses on the most used indices in LSP investigations. Specifically, this study included only vegetation indices that have been successfully tested by more than five LSP studies. These indices include but are not limited to the Normalized Difference Vegetation Index (NDVI) [41], Enhanced Vegetation Index (EVI) [42], two band Enhanced Vegetation Index (EVI2) [43], Normalized Difference Water Index (NDWI) [44], Wide Dynamic Range Vegetation index (WDRVI) [45] and recently, new spectral VIs such as the Normalized Difference Phenology Index (NDPI) [46]. Despite their successful application in LSP studies, the utility of these VIs is influenced by various factors such as sensor degradation, atmospheric impurities and snow affecting the quality of the data in the time series [18,47,48]. To address these challenges, many LSP algorithms have been developed for noise reduction, gap filling, data smoothing as well as for retrieving vegetation phenological parameters from satellite data [35,49,50]. However, there is no universal approach for estimating LSP that can be relevant in all applications. Although certain algorithms and data set combinations may produce good results in a specific land cover type, they may perform poorly in other LSP related applications. In this regard, the choice of data sets, sensors and processing methods depends on the feasibility and objectives of the LSP study.

Some of the previous LSP-related reviews focused on analyzing LSP as an indicator for global terrestrial ecosystem dynamics [51]. The review provided a synthesis of the contributions of the global LSP to the development of environmental knowledge. LSP progress and challenges in temperate and boreal forests were also reviewed recently [52]. The study uncovered an in-depth intercomparison of in situ and satellite phenological metrics in the boreal forest. A well detailed synthesis of the LSP phenological extraction methods using multispectral satellite data was recently published [21]. The review mainly uncovered the advantages and shortcomings of phenological metrics extraction methods. A sensor-specific review investigated the utility of Sentinel-2 data in phenological research, unpacking the potential and drawbacks of the data set in LSP investigations [34]. Their study discussed LSP developments using Sentinel-2 data only. However, preceding LSP reviews have two important aspects that require further attention. Firstly, the reviews focused more on croplands and forests ecosystems; a comprehensive understanding of LSP developments in rangelands has largely remained elusive. Secondly, the phenological extraction methods covered in the previous reviews focused on the developments of various phenological extractions methods, with a particular focus on their strengths and limitations in LSP studies. However, the applications and suitability of LSP phenology packages that allows the processing of satellite time series data in a single environment was not extensively covered.

To the best of our knowledge, there is no comprehensive assessment of the progress and challenges in LSP studies with a key focus on rangeland environments. Consequently, there is a need for a state-of-the-art review to establish the milestones that have been achieved in modeling LSP in rangelands ever since the availability of satellite data in the early 1970s. Therefore, this study aims to review the progress and challenges in modeling LSP in tropical rangelands. The paper initially provides a detailed account of the methodology used in selecting the literature examined in the study. Next, a detailed overview of the satellite sensor developments and associated VIs in LSP studies is provided. The paper also interrogates the commonly used and recently developed LSP data processing software packages while proposing future research directions on the remote sensing of LSP in rangeland ecosystems.

2. Literature Search and Selection of Sources

The literature search was conducted using the Scopus and Web of Science electronic scientific databases. The search terms used were ‘Land Surface Phenology’, ‘Remote Sensing’, ‘Rangelands’, and only peer reviewed LSP literature that focused on rangelands between 1972 and 2021 was considered. Specifically, LSP studies from grasslands, savannah, shrubland, woodlands land cover types were included in this study. A total number of 120 English language publications from peer reviewed journals were retained from the search process. More studies were added to the gathered literature by reviewing the publications found in the reference list of the initially retrieved sources. A total of 37 articles from the references were included, bringing the total number of sources used in this study to 157. Because the majority of LSP studies cover a wide range of land cover types, studies that covered rangelands were included in the review. The retrieved literature sources comprised articles, reviews, book chapters, conference papers and letters.

3. Satellite Sensor Developments in LSP Studies

Before the period of satellite sensors, vegetation phenology monitoring in rangeland ecosystems mainly relied on field surveys that collected phenological events of vegetation. Precisely, these surveys recorded distinct plant-specific biological changes such as the emergence of first leaves, leaf expansion, flowering and fruiting [52,53,54,55]. However, the major limitation of these field methods is the spatial extent to which the plant phenological events were collected. The plant phenological events were mostly collected at local scales by volunteers, farmers and botanists.

Table 1. Summary of widely used satellite sensors used in estimating LSP in rangeland ecosystems.

Satellite Sensor	Spatial Resolution (m)	Spectral Bands	Swath Width (km)	Acquisition Frequency	Reference
AVHRR	1100	5	2600	12 h	[58]
MODIS	250–1000	36	2330	Daily	[59]
VIIRS	375	22	3060	12 h	[60]
Landsat ETM	30	8	185	16 days	[61]
SPOT VGT	1015	4	2250	Daily	[62]
MERIS	300	15	115	3 days	[63]
Sentinel-2	10–60	13	290	5 days	[31]
PlanetScope	3–5	4	475	Daily	[64]
Himawari	500–2000	16	1000	10 min	[32]
SEVIRI	3000	12	980	15 min	[60]
Advanced Baseline Imager	500–2000	10	1000	15 min	[65]

contains a summarized chronology of polar orbiting and geostationary satellite sensors commonly used in estimating and monitoring LSP in rangeland ecosystems. The spatial, temporal and spectral characteristics of the sensors are also presented while the suitability of each of the sensors in estimating LSP in rangelands is discussed in the text. In the last decades, satellite sensors have proven to be powerful and cost-effective tools for characterizing the phenological trends of vegetation at a larger scale [21,55,56,57].

Earlier LSP studies in rangelands were pioneered by Landsat series satellite sensors with a 16 day revisit time at 30 m spatial resolution [66]. Since its official launch in 1972, the Landsat program provided a series of earth observing satellite missions with global coverage free of charge. The Landsat data archive provides rangeland managers an opportunity to assess the long-term phenological changes (multidecadal studies) and how they affect rangeland productivity at local and regional scales. Although Landsat satellites pioneered LSP investigations, only 15% (Figure 1) of the sources used in this study explored the capability of Landsat data in rangelands. This could be explained by the sensor’s low temporal resolution, which makes it ineffective in LSP applications because some plant phenological cycles change quicker than Landsat’s 16 day revisit cycle. Additionally, the presence of cloud contamination in some of the satellite images further reduces the quantity and quality of satellite data available for the characterization of phenological cycles of vegetation [35,67,68].

In the early 1980s, LSP studies shifted more towards the use of the Advanced Very High Resolution Radiometer (AVHRR) sensor. The AVHRR data became one of the most widely used data sets for LSP investigations in rangelands around the world, especially for large scale applications [58,69,70], and has been proven to be well suited for long-term phenological studies [71]. Evidence from the literature revealed that about 21% of LSP investigations in rangelands used the AVHRR data. The successful retrieval of phenological cycles of vegetation in rangeland ecosystems requires sensor platforms with a high quick return that is sufficient to capture the rapid changes of vegetation activity. Ample credit has been given to the AVHRR satellite platform for providing long term archives of 1–8 km spatial resolution of data with global coverage. In China, AVHRR data was used for the estimation of spring vegetation green-up, and the study reported that the onset of green up advanced at a rate of 0.4 to 1.9 days per decade [72]. Another study in Central Europe reported the successful application of AVHRR in multidecadal studies; the phenological patterns of vegetation revealed a general shift to an earlier start of the season (−0.54 days per year) and extended growing seasons (0.96 days per year) in Central Europe [73]. However, the applications of AVHRR data in estimating and monitoring LSP in rangelands presented several challenges that primarily originated from the instrument’s characteristics and sensor design. These characteristics include poor radiometric calibration, geolocation errors and broad spectral channels [74].

In the late 90s, the Satellite Pour l ’Observation de la Terre (SPOT) vegetation was launched with a 10 day temporal resolution. About 10% of the research studies have used the SPOT vegetation data for LSP applications in rangelands [41,75,76]. Different from other scanner sensors like AVHRR, the SPOT instrument uses the linear array system which facilitates the production of good quality imagery at a coarse spatial resolution while it maintains a significantly reduced distortion. The utility of the SPOT vegetation satellites in retrieving LSP in rangelands is, however, limited by its low temporal resolution and atmospheric impurities such as clouds and snow. Due to the 10 day repeat cycle, the SPOT satellite’s probability of capturing cloud free satellite images over an entire growing season is reduced, leading to high temporal gaps in the time series. In 2002, The European Space Agency (ESA) deployed the Medium Resolution Imaging Spectrometer (MERIS) instrument with global coverage of three days at 300 m spatial resolution. Although the instrument was mainly designed for ocean studies, 6% of studies have used the instrument in retrieving LSP in rangelands [63,77,78,79]. Unfortunately, the mission ended in 2012, limiting long term monitoring of LSP. It was never possible to conduct multi-decadal LSP investigations using the MERIS instrument.

In the early 2000s, the Moderate Resolution Imaging Spectroradiometer (MODIS) satellite was launched, collecting daily reflectance data at 250 to 1000 m spatial resolution. More than 60% (Figure 1) of the LSP studies in rangelands used MODIS data, making it the most used data set in LSP applications globally. Evidence from literature revealed that the availability of the MODIS sensor with substantively improved sensor characteristics enabled the estimation and monitoring of vegetation phenological cycles at various scales [80,81,82]. Because the MODIS instrument is now aging and its lifespan ending, the Visible Infrared Imaging Radiometer Suite (VIIRS) was launched towards the end of 2011 to enable the continuation of the MODIS data provision mission [72,83]. However, it is also imperative to note that there are technical errors linked to switching from the MODIS instrument to VIIRS such that users should consider adjustments to ensure consistent and good quality time series before phenological metrics can be extracted.

In 2015, ESA launched the Sentinel-2 satellite with a five day temporal resolution at 10–60 m spatial resolution. The applications of Sentinel-2 data in LSP have been gaining traction lately; 8% of the studies have tested the utility of the data set in rangeland applications [31,35,82]. At 10m spatial resolution, Sentinel-2 presents a huge potential for estimation of plant specific phenological cycles, which was previously difficult using early satellite sensors such as AVHRR, Landsat and MODIS. The estimation of plant specific phenological cycles, especially alien invasive species, is crucial for better management of rangelands. Nevertheless, because Sentinel-2 was launched in 2015, it lacks a global historical coverage of vegetation activity and possibly a greater impact of the sensor’s contribution in LSP monitoring will be more evident in the next few decades. The new Sentinel-3 satellite also provides a greater potential for estimating and monitoring LSP on a global scale. The 1270 km swath width of the Sentinel-3 is well suited for global applications compared to its predecessor.

Recently, Planet Labs, an aerospace company, launched PlanetScope, a constellation of more than 130 small satellites collecting multispectral images in four bands at 3 m spatial resolution daily [84]. The utility of PlanetScope data in estimating the phenology of short vegetation cycles was tested in a Kenyan rangeland [1]. Their findings highlighted that the PlanetScope data set had more cloud free observations as compared to the Sentinel-2-time series, resulting in more reliable vegetation spatial patterns as compared to Sentinel-2. The PlanetScope data has a better quick return time which is sufficient to swiftly capture developing vegetation phenological transitions such as green up and peak. Like Sentinel-2, at 3 m spatial resolution, PlanetScope data presents a more refined potential for the estimation of plant specific phenological transition dates. However, the acquisition of PlanetScope data from Sentinel Hub is not straightforward as compared to other freely available data sets such as Sentinel and Landsat that can be openly downloaded from the public data archives such as earth explorer, Google Earth Engine and ESA. To obtain the data free of charge, the PlanetScope data acquisition process requires the submission of proposals and sponsorship applications via the Network of Resources (NoR) website. Consequently, the applicability of PlanetScope data sets in rangelands monitoring remains a challenge, especially in developing countries.

In addition to the previously discussed polar orbiting satellite sensors, high temporal resolution geostationary sensors have been increasingly used for estimating LSP in rangeland ecosystems. The most used sensors include the Himawari Imager [32], Spinning Enhanced Visible and Infrared Imager (SEVIRI) [85] and Advanced Baseline Imager (ABI) [65]. The new generation of geostationary sensors have the capability of imaging the earth at 10–15 min intervals, and they have strategically positioned spectral bands that are appropriate for deriving a wide range of VIs, thus holding a huge potential in estimating the phenological cycles of vegetation in rangeland ecosystems. Because they have high frequency revisit time, geostationary sensors provide high cloud free observations, a key attribute that is crucial in tracking the phenological developments of vegetation when monitoring rangeland ecosystems in cloud contaminated regions. Despite the progress made by geostationary satellites in LSP studies, the sensors were recently launched and therefore do not have the data archives needed for the retrogressive extraction of LSP metrics. Additionally, the geostationary sensors acquire data at low spatial (3–5 km) resolution as compared to polar orbiting satellite such as MODIS, a setback that compromises the extraction LSP metrics in rangeland ecosystems at a local scale. Combined, less than 5% of studies used geostationary sensors to extract LSP metrics in rangeland ecosystems. Most of these studies were focused on specific areas of individual countries, primarily in Japan, China and the United States of America.

4. Vegetation Indices and Biophysical Variables in LSP

Table 2. Summary of widely used VIs in LSP and their formulation.

Vegetation Index	Formulation	Characteristics and Applications	Reference
NDVI	$\frac{NIR-Red}{NIR+Red}$	Large scale vegetation assessments, related to canopy structure and canopy photosynthesis	[86]
PVI	$\frac{NIR-a1\times Re{d}^{-a2}}{\sqrt{1+a_1^2}}$	Characterizes vegetation biomass and filters the effects of soil background	[87]
SAVI	$\frac{\left(1+L\right)\left(NIR-Red\right)}{\left(NIR+Red+L\right)}$	Improves NDVI sensitivity to soil background effects	[88]
NDWI	$\frac{NIR-SWIR}{NIR+SWIR}$	shows sensitivity to the changes in leaf water content	[89]
EVI	$G\frac{NIR-Red}{NIR+C_1 Red-C_{2 } Blue+L}$	Optimized to enhance sensitivity in high biomass environments	[42]
WDRVI	$\frac{a\times NIR-Red}{a\times NIR+Red}$	Enhances the dynamic range for high biomass regions	[45]
MTCI	$\frac{R_{753.75}-R_{708.75}}{R_{708.75}-R_{681.25}}$	Correlates strongly with chlorophyll content	[90]
EVI2	$G \frac{NIR-Red}{N+\left(6-\frac{7.5}{C}\right) R+1}$	Enhances the dynamic range for high biomass regions without the blue band	[43]
GRVI	$\frac{Green-Red}{Green+Red}$	Sensitive to land cover changes	[91]
NDPI	$\frac{NIR-\left(\propto \times Red+\left(1-\propto \right)\times SWIR\right)}{NIR+\left(\propto \times Red+\left(1-\propto \right)\times SWIR\right)}$	Sensitive to changes in snow cover	[46]
NDII	$\frac{NIR-SWIR}{NIR+SWIR}$	Sensitive to soil moisture storage	[92]
PPI	$-K\times In \frac{\left(NIR-Red\right)\max -\left(NIR-Red\right)}{NIR-Red)\max -\left(NIR-Red\right)soil}$	Detection of snow seasonality	[93]
PI	$\{\begin{array}{c}0, if NDVI or NDII<0\\ NDV{I}^2-NDI{I}^2\\ 0, if PI<0\end{array}$	Reduce the effects of soil and snow cover	[94]

shows a summary of the commonly used spectral VIs in LSP and their formulations. To assess vegetation cover dynamics, scientists developed spectral VIs derived from spectral data (42, 43, 92). VIs are defined as mathematical calculations of various spectral bands that are in most cases located in the visible and NIR parts of the electromagnetic spectrum (38). The time series composites of VIs computed using satellite data have provided long term global coverage archives of valuable information on vegetation activity (92–95). The NDVI pioneered the estimation and monitoring of LSP on a global scale in the early 1970s (96). Over 62% of the studies have tested the utility of NDVI in modeling LSP in rangeland ecosystems. NDVI quantifies vegetation by measuring the difference between visible red and NIR bands that are usually part of the most common multispectral sensors. The global coverage of NDVI could be useful for predicting the ecological impacts of environmental change in rangelands at various scales (97). However, it can be argued that NDVI is not effective when extracting the start and end of the season, especially in rangelands characterized by high levels of snow cover because the onset of the NDVI normally coincides with the start of the snowmelt (98). For this reason, the NDVI becomes less sensitive to slight changes in dense vegetation canopies and more sensitive to snow.

The normalized difference infrared index (NDII) was developed for the detection of vegetation water stress using the shortwave infrared reflectance, which was reported to be negatively correlated to with leaf water content because of the leaf’s absorption capacity [95]. The SWIR band is sensitive to land surface moisture; hence, many studies have used the NDII for LSP investigations in rangelands at various scales [94,96,97]. The Normalized Difference Water Index (NDWI) was proposed to be optimal in retrieving SOS metrics for areas with large snow cover [89]. The NDWI utilizes the reflectance in the NIR and Shortwave (SWIR) sections, hence it becomes ideal for discriminating green up stage from snowmelt. NDVI and NDII were combined to develop the phenology index (PI) [94]. Specifically, the PI was developed to allow the capturing of subtle changes in the SOS and EOS transition dates. The combination of two or more vegetation indices improves accuracy in LSP retrievals because the amalgamation compliments the drawbacks of each vegetation index. Although some VIs are calculated from the same spectral bands, the phenological metrics retrieved by these indices may exhibit different trends. For example, a study in Mongolia compared SOS and EOS estimated from NDVI and Simple ration (SR) which were calculated using NIR and Red spectral bands from MODIS data [98]. Their findings reported that different mathematical expressions used in these two indices would lead to the difference in vegetation phenology trends.

The Enhanced Vegetation Index (EVI) was developed as a typical satellite vegetation index for the MODIS sensor [42]. EVI was the second most used vegetation index in estimating LSP trends in rangelands. Originally, the EVI was designed to improve vegetation sensitivity in areas with high biomass and enhanced LSP estimation through de-coupling of the canopy background signal and a reduction in atmospheric influences. Following the successful applications of EVI in vegetation monitoring, many studies have used EVI time series in LSP applications [56,99,100,101]. Evidence reported from a plethora of LSP literature shows that the EVI has proven to be valuable in enhancing one-dimensionality with vegetation biophysical properties as well as reducing saturation effects normally experienced in high biomass areas, frequently experienced when using NDVI [102,103,104]. However, the blue band is a pre-requisite for the computation of EVI, and some sensors do not have the blue band, limiting the consistency of EVI across different sensors. Subsequently, a two-band combination Enhanced Vegetation Index (EVI2) was developed for sensors such as AVHRR without the blue band. The applications of EVI2 in LSP studies provided a long term EVI record that could potentially complement the NDVI record at a global scale. Several studies have reported the best similarity of the EVI2 with EVI when extracting LSP metrics. The large scale application of EVI2 reported that the phenological patterns from EVI and EVI2 data did not show any significant differences at a global level [43]. The need to provide continuity with past and present sensors has been a critical topic in LSP literature over the past few decades. Interestingly, the EVI2 came in to bridge the gap and provided continuity of LSP data across sensors.

The Wide Dynamic Range Vegetation Index (WDRVI) is a simple modification of the widely used NDVI by enhancing the dynamic range while using the same spectral bands as NDVI [45]. The WDRVI improves the robust characterization of crop physiological and phenological trends. Several other studies have used the WDRVI in LSP estimation [105,106] and have reported this index reduces saturation in high biomass regions, which is a challenge that is commonly encountered by NDVI. However, when biomass is low, NDVI remains a better choice for vegetation characterization [48]. The Green-Red Vegetation Index (GRVI) has also been an important index in LSP investigations [91,107] and has been proven to show distinct changes in vegetation even in ecosystems that show subtle changes in plant phenological appearance. The Soil-Adjusted Vegetation Index (SAVI) was mainly developed to reduce sensitivity to environmental factors such as effects of soil background on Vis [88]. Evidence from the literature has shown that the applications of SAVI in LSP are mostly suitable for areas with less vegetation cover, where the influence of soil brightness is high [108,109,110]. Similarly, the Perpendicular Vegetation Index (PVI) was also designed to decrease the sensitivity to soil reflectance; it has a higher signal-to-noise ratio compared to NDVI [87].

The MERIS Terrestrial Chlorophyll Index (MTCI) is a sensor-specific spectral vegetation index that has been used by numerous studies in LSP investigations [63,79,110,111]. Evidence gathered from the studies shows that the MTCI is sensitive to chlorophyll changes, making it appropriate for phenology investigations. However, because the MERIS sensor ceased operation, no further developments or studies can be continued using the MTCI, a major disadvantage of sensor specific spectral indices. With the launch of more satellite sensors such as Sentinel-2, the number of spectral bands has increased, leading to the development of new VIs that can be widely utilized in LSP applications at various scales. New VIs are gaining traction in LSP studies. The Normalized Difference Phenology Index (NDPI) is a spectral index that is less affected by snow when extracting green up date [48]. The new NDPI significantly minimize the influence of snowmelt on retrieving the phenological metrics ecosystems. The Plant Phenology Index (PPI) was designed to untangle vegetation start of season and snow seasonality [93]. The PPI significantly improves the retrieval of LSP metrics especially in regions with seasonal snow cover.

While VIs remains the most used vegetation indicators, some vegetation physiological parameters have also been used to extract LSP in rangeland ecosystems. The most commonly used biophysical variables include the leaf area index (LAI) [112,113,114], Fraction of Absorbed Photosynthetically Active Radiation (FAPAR) [115,116], Solar-Induced Chlorophyll Fluorescence (SIF) [117,118] and Vegetation Optical Depth (VOD) [22]. Evidence from the literature shows that the parameters provide detailed information about the biophysical characteristics of vegetation and how they respond to climatic and environmental changes [112]. The LAI biophysical variable has been widely used in LSP applications because of the availability of numerous remote sensing products and its clear signal of the physical and biological processes that are related to vegetation phenological cycles at various scales [119,120].

5. LSP Software Packages for Data Processing

The spectral reflectance from vegetation can be disturbed by a variety of atmospheric impurities such as ground and atmospheric conditions [121], changes in the satellite sensor’s illumination patterns and viewing angle [122], causing inaccurate trends in the time series data. Figure 2 shows a detailed flowchart of LSP procedures. The widely used phenological metrics are shown in a schematic diagram under the phenological metrics retrieval section. The successful estimation of LSP requires three main steps in processing spectral VIs time series data. These stages include (1) preprocessing of VIs time series data by detecting and removing outliers, (2) data smoothing and gap filling, and (3) extraction of LSP metrics. The performance and limitations of LSP data smoothing methods have been extensively reviewed [21,123]. Therefore, the developments and drawbacks of data smoothing techniques and phenological extraction methods will not be included in this section. The focus will be on reviewing the widely used LSP software packages that are designed to provide functions that perform data smoothing, gap filling and LSP metrics retrieval in a single environment. The drive for developing such all-in-one LSP software packages was to avoid errors related to moving large amounts of time series data from one algorithm to another for further processing.

The Harmonic Analyses of NDVI Time-Series (HANTS) algorithm was proposed for the extraction of the characteristics of the vegetation dynamics. The HANTS algorithm has been extensively used in rangeland LSP investigations successfully [111,124,125,126]. The amplitude and harmonic components of the algorithm makes HANTS attractive for LSP studies [125]. However, evidence from the literature shows that HANTS limits the ability of users to specify phenological parameters and it also preserves year to year variations in the time series data as it filters atmospheric impurities such as clouds. The TIMESAT program provides tools for analyzing time series data for various LSP applications [127]. The TIMESAT program has been extensively utilized by researchers to estimate and monitor LSP in rangelands [86,128,129,130]. The TIMESAT program offers three data filtering functions to remove noise from the time series data. These are asymmetric Gaussian (AG), double logistic (DL), and adaptive Savitzky-Golay (SG) filtering. One of the major advantages of using TIMESAT is that it enables the flexible tuning of thresholds and parameter settings compared to other algorithms such as HANTS. As a prerequisite, the TIMESAT program requires gap free time series data for successful characterization of the phenological cycles of vegetation cycles. The gap free pre-requisite limits the applications of the software package especially in circumstances where the data is of poor quality due to persistent cloud contamination. Consequently, TIMESAT was revised to produce temporally and spatially continuous time series data [131]. The revised TIMESAT was later named enhanced TIMESAT in a study that estimated vegetation phenology metrics from MODIS data [130]. The enhanced TIMESAT significantly improved the time series data quality by replacing the upper envelope function in the original TIMESAT software with the tuned data quality weightings whereby higher quality retrievals were assigned more influence than low quality retrievals.

TimeStats is a software package used for retrieving temporal patterns of vegetation activity from satellite data [132]. Unlike other packages such as TIMESAT that use various least-squares approaches for preprocessing time series data, TimeStats provides robust statistical methods such as cyclic, transient and stochastic for preprocessing of time series data using a pixel-by-pixel approach. Phenosat was developed as a semi-automated tool for the temporal analysis from satellite data [133]. The Phenosat tool can capture double growth season and is flexible in the selection of data interval process, a limitation that was encountered using other software such as TimeStats and TIMESAT. Although TIMESAT, HANTS and Phenosat are freely available for use, they require Matlab software, which is a nonfree environment to execute the analysis and therefore limits their applications in resource constrained regions. An open source QPheno-Metrics software package was proposed for extraction of vegetation phenological metrics, which is executed in QGIS, a Geographic Information Systems (GIS) software that is freely available [49].

Recently, the utility of the multisource land surface phenology (MS-LSP) algorithm that combined data from fine spatial resolution Sentinel-2A and -2B (10 m) and medium spatial resolution Landsat 8 (30 m) was tested [35]. The launching of the MS-LSP algorithm is a significant achievement in LSP investigations because the temporal frequency of Landsat 8 was not sufficient for the estimation of vegetation phenological changes. The current study shows that the successful estimation and monitoring of LSP in rangeland ecosystems relies heavily on the availability of robust algorithms that are capable of processing vegetation time series while minimizing atmospheric noise and sensor-related errors. In developing countries, the availability of these software packages free-of-charge plays a crucial role in data modeling for the management of rangelands. Concerning the choice of LSP algorithm to use, it is important to note that each method has its strengths and weaknesses. Therefore, the suitability of each smoothing algorithm largely depends upon the objectives of that study and the targeted land cover types. The level of expertise demonstrated by the user in tuning parameters and settings during the data processing and phenology metrics extraction process influence the results.

6. LSP Metrics Validation

To evaluate the accuracy of remotely sensed time-series data in estimating the phenological cycles of vegetation in rangeland ecosystems, there is a need for comparisons between satellite-based phenological metrics and what can be observed or measured on the ground [22,134,135]. The linking of ground-based phenological observations with satellite-based LSP retrievals provides a huge potential for tracking the response of vegetation to climatic changes in rangelands [136]. However, large samples of ground phenological observations or measurements are required to adequately validate LSP retrievals, especially in cases where studies are conducted at a large scale. The widely used methodological approaches for LSP validation include the periodic recording and documentation of distinct plant specific biological changes such as the emergence of first leaves, leaf expansion, flowering and fruiting at a local scale. The phenological stages are usually recorded by farmers, botanists, naturalists, volunteers as well as phenological observation networks such as the Nature Canada PlantWatch [137] and the USA National Network (USA-NPN) [138]. For hundreds of years, botanists and naturalists have been collecting a diversity of plant specimens in the world’s herbaria [139,140,141]. Eventually, scientists and researchers recognized the potential use of herbarium specimens in detecting and modeling the phenological changes of plants in response to environmental and climatic changes [142,143,144]. However, the collection of herbarium specimens data is labor intensive and time consuming, and this may lead to inconsistencies in the data collection protocols and methodologies.

Some researchers and scientists defined the first appearance of green leaves as the green up date [135,145,146], which resembles the start of the season in LSP. However, due to scale mismatch and data uncertainties, the phenological phases observed and recorded on the ground may sometimes be inconsistent with satellite-based LSP estimations [40]. Furthermore, field-based visual surveys have been reported to be costly and time-consuming [106], making it difficult to effectively use them in validating large scale LSP retrievals. Digital repeat photography widely known in remote sensing literature as Phenological Cameras (PhenoCams) have also been extensively used in validating LSP retrievals. For example, a study in Northern Japan used PhenoCam data from the Phenological Eyes Network to evaluate the LSP metrics from the Himawari Satellite Imager and MODIS [147]. Their study highlighted that the Root Mean Square Difference (RMSD) between LSP metrics from MODIS and PhenoCam data was high in spring and fall (14–29 days). Similarly, another study in the United States reported that LSP metrics estimated using VIIRS and MODIS showed an agreement with SOS and EOS metrics estimated from PhenoCam data with RMSDs ranging from 10.1 to 21 days [83].

However, despite the successful application of PhenoCams in retrieving phenological metrics that are comparable to satellite-based retrievals, several other studies reported low or no correlation between the two data sources, especially in rangelands [148,149]. Generally, findings from the literature show that the inconsistencies that may arise between PhenoCams and satellite-based LSP metrics are largely linked to discrepancies caused by differences in the scale of observation. In some cases, the phenological changes of vegetation in small areas monitored by PhenoCams cannot accurately be representative of large-scale changes observed by satellite images. Other differences could emanate from the fact that PhenoCams do not acquire images in the NIR wavelengths [150], they only rely on information from the visible bands. Unlike satellite sensors which have global coverage, PhenoCams coverage is limited. The majority of PhenoCam networks have specific regions they cover, for instance, Phenological Eyes Network mainly covers Japan, USA, China, Malaysia, and the United Kingdom [147], while the PhenoCam network (http://phenocam.sr.unh.edu: Accessed on 14 January 2021: 18:09) mostly focuses on the terrestrial ecosystems of North America [151]. Consequently, there is little or no PhenoCam coverage in rangelands that are geographically located in poor regions, making it difficult to validate LSP retrievals in those under-resourced areas.

Unmanned aerial vehicles (UAVs), also referred to as drones in the literature, have been an important data source for validating satellite-based phenological metrics in many rangeland ecosystems [52,152]. Carbon dioxide flux observations from FLUXNET stations have also been used to validate LSP retrievals in rangelands. For instance, a study in Canada compared LSP retrievals from MODIS and SPOT with global FLUXNET-derived phenological data [41]. Their study reported that the satellite-based LSP estimations had an overall low correlation (R² < 0.30) with the phenological timings obtained from the flux observations. On the contrary, another study highlighted that carbon flux phenology estimations were highly comparable to satellite-based LSP metrics, with R² values ranging from 0.43–0.78 amongst the four biomes used in their study [149]. However, the use of carbon flux observations to estimate the phenological cycles of vegetation is more suitable in biomes with distinct and detectable seasonal cycles because it can be arduous to detect the start of carbon uptake in high biomass environments [148]. As widely reported in the literature, another limitation associated with the use of carbon flux data in validating LSP retrievals is the limited coverage of eddy covariance flux data collection sites in many biomes [149,153,154], thus the estimation of vegetation phenological cycles using carbon flux data remains challenging at a large scale. Another method of validating LSP retrievals involves the use of bioclimatic models which use precipitation and temperature data to track the phenological changes of vegetation [155,156]. However, several studies have reported low correlations between satellite-derived phenology and bioclimatic-based phenology retrievals [157,158]. Evidence from the literature consulted in this review shows that the majority of the LSP validation methods lack detailed spatial and temporal ground phenological measurements and events which include species level field observations.

7. Challenges and Future Directions in Rangeland LSP

LSP research has been hobbled by inconsistencies between remote sensing retrievals and vegetation phenological events recorded on the ground [83,159,160]. Evidence from the literature has shown that the sources of disagreements in phenological metrics amongst different vegetation species may arise due to VIs used [98,160], satellite sensor characteristics [24,72], atmospheric conditions [161] and algorithms used to smooth and estimate the LSP phenological metrics [18,24,162]. The variations introduced by atmospheric impurities, absorption, and scattering can considerably affect the precision of users when interpreting remote sensing images, especially for the detection of vegetation dynamics at a landscape scale [163]. Additionally, the impacts of bandpass and other sensor characteristics on the behavior of the VIs causes errors when retrieving LSP metrics. Another common challenge in LSP investigations is mixed pixels. A pixel in the VIs time series may contain an unknown composition of vegetation species and may result in mixed signals [164] because vegetation species vary in their sensitivity to climatic fluctuations and changes.

Although it does exist, future LSP research should probably consider investing more in developing algorithms that fuse moderate and high spatial resolution sensors to improve LSP metrics retrieval at the species level. The extraction of phenological metrics in rangelands such as alpine grasslands remains a challenge due to subtle seasonal variation in VIs time series [109,165,166]. To address saturation problems encountered when deriving phenology metrics, there is a need for further research on the reconstruction and regeneration of VIs that will adequately reduce the problem of saturation in high biomass regions. Specifically, the WDRVI was reported to be effective in dealing with saturation, a common problem encountered using NDVI [45]. Ground-based phenological observations provide reliable and accurate information on individual plant species, but LSP observes changes at a large scale; hence, the use of species-specific phenological phases observed on the ground to validate large-scale satellite LSP retrievals becomes problematic. The current study suggests that international phenology research networks such as International Long Term Ecological Research Network (ILTER) have the potential to facilitate the regulation and standardization of phenology research protocols.

The compatibility of formats when linking LSP algorithms and remote sensing data processing software packages is a challenge. For instance, phenology packages such as TIMESAT [167] do not perform the preprocessing of satellite data and the computation of vegetation indices. Consequently, data conversion procedures may lead to the loss or deformation of valuable information in the time series. The development of LSP algorithms that can preprocess remotely sensed data and handle phenological analysis in the same environment will go a long way in addressing some of the technical errors encountered in LSP retrieval. Specifically, deep learning presents an opportunity to further develop robust software packages [168] with the abundance of data processing tools and techniques that can be used to better characterize the phenological cycles of vegetation in rangelands. The validation of LSP metrics remains a challenge in many of the biomes covering the African continent.

To improve the LSP estimation and monitoring, there is a need for the establishment of more phenology networks in Africa. While a plethora of LSP studies has focused on the investigation of croplands and forests, the study of the phenological cycles of alien invasive species in rangelands has largely lagged. The invasive species LSP modeling will enable understanding of the biological structure and timing of alien invasive phenology and lead to improved rangeland management using this knowledge to choose suitable treatment methods in zones infested by the alien invasive plants. The applications of medium spatial resolution sensors such as Landsat in LSP studies have been largely limited by the sensor’s temporal resolution because some plant phenological appearance changes quicker than Landsat’s 16-day temporal resolution. On the other hand, the high temporal resolution sensors such as MODIS lack spatial detail that can accurately track the phenological events of vegetation at the subspecies level. To tackle these challenges, the constellation of Planet sensors (Planet Scope) holds a huge potential in LSP applications. PlanetScope data have an appropriate quick return time to capture swiftly developing phenological transitions, such as green-up.

The effective management of rangelands requires timely, accurate and reliable information about vegetation activities and how they change over time due to changes in climate or environmental conditions [11]. Based on the LSP developments discussed in this study, there is a huge potential of LSP data for future use in a wide range of applications for better management of rangeland ecosystems. User groups such as rangeland managers, ecologists and government agencies may use published LSP data sets for various objectives such as designing better management systems and planning restoration activities after ecosystem disturbances. A typical example of the relevance of LSP data is the forecasting of grass and invasive species reproduction based on phenological models that incorporates satellite data, climate data and phenological networks observation networks databases.

8. Conclusions

The current study has reviewed the literature on the progress of remote sensing in estimating land surface phenology in rangelands. Empirical evidence has shown that remote sensing offers invaluable data sources in the generation of phenological trends of vegetation at regional and global scales, a task that was previously impossible using ground-based phenological observations. The application of remote sensing in LSP studies was pioneered by early satellites such as Landsat, AVHRR and MODIS. However, evidence from the literature revealed that the applications of these data sets have been largely limited by the sensor’s low spatial resolution, poor radiometric calibration and geolocation errors. Consequently, the launching of high spatial resolution sensors such as Sentinel-2 has significantly improved rangeland LSP investigations. Findings from the literature revealed that NDVI pioneered LSP investigations in the early 1970s, and many other indices were subsequently developed to address the shortcomings of NDVI. Although milestones have been achieved in the applications of VIs in the retrieval of LSP, the modeling of phenology from remote sensing remains a challenge because it is difficult to develop VIs models that can be used efficiently in all environments. Owing to the unique traits of many Vis and models, appropriate indices should be used to characterize vegetation phenological cycles at various growth stages and estimate phenology trends in different biomes. The successful retrieval of LSP metrics depends on the availability of robust algorithms that are capable of processing vegetation time series while minimizing atmospheric noise and sensor-related errors. In developing countries, the availability of these software packages free-of-charge plays a crucial role in data modeling for the management of rangelands. Currently, there is no agreement in the LSP community regarding the best methods and algorithms for retrieving LSP metrics. In the future, LSP researchers should invest more in the development of LSP algorithms that can preprocess remotely sensed data and handle phenological analysis in the same environment as well as facilitating the regulation and standardization of phenology research protocols for better management of rangeland ecosystems.