Leveraging Sentinel-2 and Sentinel-1 imagery for retrieving key wheat phenology stages in a dryland agricultural setting

Monitoring crop phenology is crucial for agricultural management, providing the basis for timely decisions related to irrigation and fertilization, while offering insights into crop condition. This is particularly important in drylands, where agriculture is constrained by water scarcity and harsh environmental conditions, making phenology monitoring essential for optimizing resources and sustaining productivity. Despite its importance, traditional ground-based phenology monitoring methods are limited in scalability and consistency, underscoring the need for satellite-based alternatives. In this study, we investigated the feasibility of using time-series data from Sentinel-2 and Sentinel-1 imagery, each analyzed separately and combined with curve-based rules (i.e., threshold- or inflection-point criteria applied to temporal profiles), to retrieve five phenology stages of germination date, flowering date, senescence date, harvesting date, and cycle length for a large number of wheat fields in the Al-Jouf region of Saudi Arabia, a hot desert environment receiving less than 100 mm of annual precipitation that typifies dryland agriculture. To achieve this, we used 85 cloud-free Sentinel-2 and 94 Sentinel-1 images, along with in situ phenology observations. From Sentinel-2 imagery, we derived the normalized difference vegetation index (NDVI), green-red vegetation index (GRVI), and plant senescence reflectance index (PSRI), while volume backscattering (VV-VH), i.e., difference between vertical transmit/vertical receive (VV) and vertical transmit/horizontal receive (VH), was derived from Sentinel-1 imagery. Field-level time-series of each Sentinel-2 and Sentinel-1 metric were extracted and fitted with Gaussian functions. Next, the five phenology stages were retrieved by applying stage-specific curve-based rules to Gaussian-fitted field-level time-series of each metric. The accuracy of retrievals was assessed against in situ observations. Results showed that curve-based rules applied to NDVI, GRVI, and PSRI provided relative root mean square error (rRMSE) values of 2.64%–14.17% and bias magnitudes of 7.34–20.60 days across all stages without requiring rule refinement. Rules applied to VV-VH produced rRMSE values of 4.60%–26.22% and bias magnitudes of 7.11–43.73 days. Targeted refinement of flowering and harvesting rules reduced these errors. Overall, this work demonstrates the effectiveness of applying curve-based rules to each of Sentinel-2 and Sentinel-1 data separately in the dryland setting examined here, providing a foundation for evaluation across the broader range of drylands.

Wheat is a staple crop that plays a vital role in global food security (Meng et al., 2009). Its importance is even greater in drylands, where agricultural productivity is constrained by climate variability and limited water availability (Naorem and Machiwal, 2023). In dryland environments, there is a pressing need for analytical insights from continuous monitoring and management strategies to optimize resource use and enhance agricultural resilience and productivity. Monitoring and understanding crop phenology, i.e., the timing of key growth stages, may offer a promising pathway to support effective agricultural management and production practices (Zhou et al., 2023). Accurate crop phenology information enables timely agricultural management decisions to be implemented, including the application of irrigation, fertilization, and pest control, tailored to specific crop growth stages (Pimstein et al., 2009). Information on crop phenology stages also serves as a valuable indicator of crop condition and productivity (Meng et al., 2009), provides insight into climate variability responses (Beaubien and Freeland, 2000; Begue et al., 2014), and supports crop growth simulation models that assess crop development under varying environmental and management scenarios (Clevers et al., 1994; Sun, 2000).

Despite the importance of crop phenology monitoring, traditional ground-based methods used for this purpose face significant limitations in terms of scalability, consistency, and objectivity. Common approaches such as visual scoring (Bogard et al., 2011), soil plant analysis development (SPAD) meter readings (Xie et al., 2016), and proximal sensing measurements (El Hajj et al., 2024) are often subjective, labor-intensive, and impractical for large-scale applications (Christopher et al., 2014). These methods typically provide observations from limited areas or a small number of plants, which may not capture the spatial and temporal variability in crop phenology across diverse landscapes. Such variability is driven by differences in weather conditions, management practices, planting dates, crop varieties, and local environmental factors (Estrella et al., 2009; Sakamoto et al., 2010). Accurately capturing this heterogeneity at large spatial scales is challenging with ground-based tools alone (Anderegg et al., 2020), underscoring the need for spatially explicit and temporally consistent phenology retrieval methods. In this context, satellite-based observations offer a scalable and consistent solution for capturing crop phenology dynamics across large agricultural areas (Liang et al., 2011; Zhou et al., 2023).

Satellite imagery provides spatially continuous and temporally frequent data at regional to global scales that can be used to support various information services for monitoring agricultural and environmental processes (Jensen, 2009). The European Space Agency launched the Sentinel missions under the Global Monitoring for Environment and Security (GMES) program (Attema et al., 2009; Spoto et al., 2012). Sentinel-1 is a C-band synthetic aperture radar (SAR) satellite constellation that offers day-and-night imaging collection that is unaffected by the majority of atmospheric conditions, making it highly suitable for use in persistently cloudy regions or in urgent monitoring scenarios (Attema et al., 2009). Sentinel-2, on the other hand, provides high-resolution optical observations through its multi-spectral instrument, which captures 13 spectral bands covering the visible, near-infrared, and shortwave infrared regions (Spoto et al., 2012). The key distinction between the two missions lies in their sensing technologies with Sentinel-1 being an active radar system, emitting its own signal, operating in all weather, and lighting conditions, while Sentinel-2 is a passive optical system that relies on reflected sunlight and is constrained by cloud cover and daylight availability.

Numerous studies have investigated wheat phenology retrieval, though many have focused on limited aspects or adopted methods with narrow applicability. For instance, Mercier et al. (2020), El Imanni et al. (2022), Lobert et al. (2023), Katal et al. (2022), and Chuine and Régnière (2017), have used satellite observations to detect multiple wheat phenology stages. However, many of these approaches require extensive calibration, large training datasets, or high computational resources, which are conditions that are difficult to meet in data-scarce or operational contexts. In contrast, Pan et al. (2015), Guo et al. (2016), Sun et al. (2024), and Lu et al. (2014) explored simpler retrieval methods, but these were generally limited to a single remote sensing data source (either optical or radar) and often focused only on detecting wheat germination and harvesting dates. This narrow focus overlooks intermediate wheat growth stages, such as flowering and senescence, which are crucial for in-season agricultural decision-making. Furthermore, relying on a single sensor limits the ability to assess the relative strengths of optical and radar data in capturing different aspects of wheat growth dynamics. Importantly, the aforementioned studies were conducted in non-dryland settings, limiting their relevance to dryland regions, where several unique challenges create a need for satellite-based wheat phenology retrieval methods tailored to dryland conditions. Unlike temperate regions, phenology timing and retrieval via remote sensing in drylands is governed by field-specific irrigation calendars, management decisions, and high inter-annual climatic variability, compounded by complex soil-vegetation spectral mixing and limited in situ data availability. These factors depart from the assumption of rainfall as the dominant phenology driver, which underlies most phenology retrieval approaches developed for rainfed regions. This distinction has direct consequences for temporal-series-based methods for retrieving phenology timing. For instance, frequent irrigation and high climatic variability can introduce irregularities into the temporal backscatter profile due to soil moisture fluctuations and water droplet on the canopy, an effect not frequently observed in rainfed areas, where backscatter varies more gradually in line with sparse rainfall. A comparable effect occurs in vegetation indices, where irrigation-driven fluctuations in canopy moisture and greenness can produce peaks from the smooth single-peaked seasonal curve typically used to define transitions in phenology in rainfed systems, complicating the direct transfer of methods calibrated for rainfed regions to irrigated dryland settings. Lastly, complex soil-vegetation spectral mixing over sparse canopies and bright soils weakens the vegetation signal in optical indices, with elevated surface scattering similarly affecting radar observations, particularly during the low-biomass periods surrounding germination and harvest (Smith et al., 2019), and limited in situ data availability constrains the use of methods requiring extensive calibration or large training datasets. In general, studies implementing operationally scalable methods for retrieving multiple wheat phenology stages remain scarce. This is especially true for approaches using curve-based rules, defined here as threshold- or inflection-point criteria applied to fitted temporal vegetation index or backscatter profiles, which require minimal to no calibration and offer a straightforward path to broad applicability.

The objective of this study is to assess the feasibility of using time-series data of optical vegetation indices from Sentinel-2 imagery and volume backscattering from Sentinel-1 imagery, combined with curve-based rules, for retrieving five key phenology stages, including germination date, flowering date, senescence date, harvesting date, and cycle length, from each satellite separately across a large number of wheat fields. To the best of our knowledge, this study is among the few to develop a rule-based retrieval approach specifically tailored for data-scarce regions. The developed approach is applied in a dryland agricultural setting where investigations on wheat phenology retrieval remain limited, and is rigorously evaluated against an extensive in situ dataset spanning two growing seasons. Furthermore, optical and radar data are used independently, enabling a direct comparison of the respective capabilities of Sentinel-2 and Sentinel-1 in capturing different aspects of wheat growth dynamics. By relying on publicly available satellite data, this study offers a practical, scalable, and resource-efficient solution for wheat phenology monitoring in dryland agriculture, reducing dependence on labor-intensive and spatially limited ground-based observations while providing valuable insights for agricultural management and productivity optimization in data-scarce regions.

Study site description

The study site is a large commercial farm located in the Al-Jouf region of northern Saudi Arabia (see Figure 1a). The climate of Al-Jouf is classified as hot desert (BWh) according to the Köppen-Geiger classification (Peel et al., 2007), with annual precipitation below 100 mm (El Kenawy and McCabe, 2015) and mean air temperatures ranging from 10 °C in winter to 40 °C in summer (Elfarkh et al., 2023), representing a dryland environment typical of many regions in the Middle East. Despite these harsh conditions, Al-Jouf is among Saudi Arabia’s most important agricultural regions, with an irrigated area that expanded rapidly from almost non-existent in the 1980s to more than 2,500 km 2 by 2015 (Lopez et al., 2020). The selected farm itself covers 30,000 ha and consists of 300 center-pivot irrigated fields (see Figure 1b). Wheat is the dominant winter crop grown on the farm and is typically sown between early November and late December each year. It is generally harvested between early May and early June of the following year. Although the farm operates under conventional center-pivot irrigation practices, spatial and temporal heterogeneity in wheat phenology at the farm scale was anticipated due to differences in field-level management across the center-pivot fields. Such heterogeneity, combined with the region’s climatic extremes, made the farm an ideal testbed for evaluating satellite-based wheat phenology retrieval approaches tailored to dryland agricultural systems. According to farm-reported soil analyses, the soil texture is sandy loam (76% sand, 14% silt, 10% clay) with a wilting point of 7.5 vol%, a field capacity of 14.8 vol%, and a saturation point of 43.5 vol%.

Data collection and description

Sentinel-2 data

The Sentinel-2 mission consisted of two operational satellites (Sentinel-2A and Sentinel-2B) during the study period (2019–2020 and 2020–2021 growing seasons), equipped with a multispectral instrument providing imagery at high spatial and temporal resolutions across 13 spectral bands (Spoto et al., 2012). This study used 85 cloud-free Sentinel-2 surface reflectance (Level-2A) images, acquired over the study site for the two growing seasons with a nominal revisit time of 5 days. The used images span a period from the beginning of October to the end of July for each growing season (see Figure 2), covering the full growth cycle of each active wheat field, from pre-sowing to post-harvest.

Three vegetation indices, including the NDVI, plant senescence reflectance index (PSRI), and green-red vegetation index (GRVI) were derived from the Sentinel-2 imagery. NDVI (Equation 1) is widely used to monitor vegetation greenness and was calculated from the near-infrared (band 08) and red (band 04) bands (Rouse, 1973). PSRI (Equation 2) is sensitive to pigment degradation during senescence and was derived from the red (band 04), blue (band 02), and red-edge (band 06) bands (Merzlyak et al., 1999). GRVI (Equation 3) captures early green-up, autumn coloration, and subtle vegetation changes or ecosystem disturbances, and was computed from the red (band 04) and green (band 03) bands (Motohka et al., 2010). All of these spectral bands have a spatial resolution of 10 m, except band 06, which has a resolution of 20 m. Band 06 was therefore resampled to 10 m using the default resampling applied by Google Earth Engine during Sentinel-2 image export, allowing all three vegetation indices (NDVI, PSRI, and GRVI) to be computed at a 10 m spatial resolution. These vegetation indices were selected to retrieve wheat phenology stages due to their complementary sensitivities to different stages of wheat (and general crop) development and their demonstrated effectiveness in wheat phenology retrieval, as reported in previous studies (Pan et al., 2015; Anderegg et al., 2020; El Imanni et al., 2022; Zhou et al., 2023; Sun et al., 2024). where B 08 is the near-infrared band; B 04 is the red band; B 02 is the blue band; B 06 is the red-edge band; and B 03 is the green band.

Sentinel-1 data

The Sentinel-1 constellation comprised two operational polar-orbiting satellites (Sentinel-1A and Sentinel-1B) during the study period, providing C-band SAR data under consistent acquisition conditions (same orbit, pass direction, and incidence angle). The Sentinel-1 observation system operates in interferometric wide swath mode as the default configuration over land, producing dual-polarized imagery, vertical transmit/vertical receive (VV) and vertical transmit/horizontal receive (VH) (Attema et al., 2009). For this study, 94 Sentinel-1 images were acquired over the study site from the beginning of October until the end of July for both the 2019–2020 and 2020–2021 growing seasons with a nominal revisit time of 6 days (see Figure 3). All used images were from relative orbit 116 with incidence angles ranging from 30° to 32.5°, and were provided in Level-1 ground range detected (GRD) format. The Sentinel-1 GRD images were provided with a spatial resolution of 10 m and were pre-processed using Google Earth Engine, following a standard pipeline consisting of: (1) application of the orbit file; (2) border noise removal; (3) thermal noise removal; (4) radiometric calibration; (5) terrain correction (orthorectification); and (6) conversion to backscatter coefficient in decibels, as detailed in (Filipponi, 2019; Tufail et al., 2022).

In this study, the volume backscattering was estimated using the difference between VV and VH polarizations (VV-VH), which enhances sensitivity to vegetation volume by suppressing surface backscattering contributions from soil moisture, texture, and roughness. Since surface scattering tends to affect both VV and VH in a quasi-similar sensitivity (El Hajj et al., 2019), subtracting VH from VV helps isolate the volume scattering signal primarily originating from the wheat (or general crop) canopy. This suppression of surface-scattering contributions is particularly relevant in the present study area, where bright sandy loam soil may otherwise contribute substantially to the backscatter signal. This differential approach results in a cleaner vegetation signal that is more responsive to structural changes, such as biomass accumulation and development, making VV-VH a suitable indicator for monitoring wheat (or general crop) phenology stages (Vreugdenhil et al., 2020; Ling et al., 2022).

In situ data

In situ data were obtained from the farm for 68 and 73 active wheat fields during the 2019–2020 and 2020–2021 growing seasons, respectively, to serve as ground truth for validating the retrieved wheat phenology stages. The collected in situ data included germination date, flowering date, senescence date, harvesting date, and cycle length. An exception applied to flowering and senescence dates, for which the collected in situ records related to these two stages were available only for 73 active wheat fields during the 2020–2021 growing season. All collected data were observed and recorded by farm staff as part of a standardized phenology monitoring system that supported both tracking wheat development and documenting key management activities such as fungicide applications and nutrient adjustments. Since wheat phenology transitions do not necessarily occur simultaneously across all plants within a field, dates were assigned through visual field assessments, when a wheat field was considered to have reached more than 50% of the plants exhibiting the corresponding phenology characteristics. This criterion aligns with common field phenology observation practices, such as those described in (Koch et al., 2007). All phenology dates were subsequently converted to day of year (DOY) format (see Table 1). According to farm management, a single wheat cultivar (i.e., Yecora Rojo) was cultivated across all fields during both growing seasons.

Summary statistics describing the variability of the collected in situ data are presented in Table 1. Germination dates showed the greatest variability across both growing seasons, likely due to staggered sowing schedules across the large number of fields, as germination typically occurs shortly after sowing and is therefore strongly influenced by sowing timing. In contrast, harvesting dates over both seasons, as well as flowering and senescence dates, exhibited narrower interquartile ranges, indicating synchronized wheat development across fields and reflecting broadly comparable agronomic and environmental conditions (see Table 1). Regarding cycle lengths, they mirrored the variability observed in germination dates.

Retrieval of wheat phenology stages from Sentinel-2 data

To retrieve the key wheat phenology stages (germination date, flowering date, senescence date, harvesting date, and cycle length) using Sentinel-2 imagery, we analyzed time-series of the three vegetation indices: NDVI, PSRI, and GRVI. NDVI was employed to derive the germination date, harvesting date, and cycle length, consistent with (Pan et al., 2015; El Imanni et al., 2022; Sun et al., 2024), who demonstrated NDVI’s effectiveness in identifying such stages. PSRI was used for identifying the senescence date, following (Anderegg et al., 2020), who showed that PSRI consistently outperformed NDVI and all other tested vegetation indices (50 in total), accurately tracking visually observed senescence dynamics across diverse wheat genotypes and environments. GRVI was selected for detecting the flowering date, as it has been shown to reliably capture early reproductive development in wheat with values peaking at the heading stage (Zhou et al., 2023), which directly precedes the flowering stage in wheat (Miller, 1992).

NDVI, PSRI, and GRVI were produced for each Sentinel-2 image using the equations outlined in Section 2.2.1. Active wheat fields were delineated and labeled in the QGIS software (v3.22) (Moyroud and Portet, 2018) for both growing seasons. The delineated field boundaries were then buffered inward by 10 m to exclude edge pixels and reduce mixed-pixel effects at the field margins. For each image, the vegetation index values of all pixels within each delineated wheat field were averaged to construct field-level time-series profiles. Subsequently, further processing steps were applied to the index time-series of each field in order to retrieve wheat phenology stages. First, each index time-series was normalized using the min-max method to rescale values between 0 and 1, according to the following equation (Vafaei et al., 2022): where N t is the normalized vegetation index value at time t; V t is the original vegetation index value at time t; and V min and V max are the minimum and maximum index values across the time series, respectively. Next, the normalized time-series data were smoothed using the Savitzky-Golay filtering technique, which reduces noise while preserving the overall shape of the curve. The filtering process was applied using the following equation (Savitzky and Golay, 1964): where Y j * is the smoothed vegetation index value at position j; Y j + i is the original vegetation index value at offset i; C i is the Savizky-Golay coefficient; and N is the number of convoluting integers which is equal to the smoothing window size 2 m + 1 (Niedźwiecki et al., 2021). In this study, a smoothing window size of 7 and a second-order polynomial were used for the Savitzky-Golay smoothing. Finally, a Gaussian function was fitted to the smoothed index time-series curves. This function was selected for its suitability in modeling the bell-shaped profiles typical of NDVI, PSRI, and GRVI seasonal trajectories. The Gaussian fitting was performed using the following equation (Luh, 2012): where f t is the fitted vegetation index value at time t; A is the amplitude (the peak of the curve); μ is the mean (the time at which the peak occurs); and σ is the standard deviation (which determines the spread of the curve over time).

Next, we employed insights from previous wheat phenology studies (Large, 1954; Miller, 1992; Christopher et al., 2014; Pan et al., 2015; Anderegg et al., 2020; Alabama Cooperative Extension System, 2022; El Imanni et al., 2022; Zhou et al., 2023; Sun et al., 2024) to design curve-based rules (i.e., threshold- or inflection-point criteria applied to fitted temporal profiles) for retrieving phenology stages from the fitted vegetation index curves. These studies employed processed time-series of NDVI, PSRI, and GRVI with bell-shaped seasonal trajectories similar to those obtained herein. Such similarity makes the transfer of established insights both feasible and reliable. The designed rules specified positions along the x- and y-axes of the fitted curves using the Gaussian parameters (A, μ and σ; see Equation 6) and were applied individually to derive phenology stages for each field. Although originally developed based on insights from studies that were conducted under site-specific conditions, the reliance of these rules on relative curve positions, together with the consistency in curve shapes, supports their adaptation to the current study context. This approach enabled the operationalization of established insights and allowed us to assess their robustness in capturing wheat phenology stages across a large number of fields. The following subsections detail how each phenology stage was retrieved from Sentinel-2 data.

Germination date retrieval

Building on the work of Pan et al. (2015) and Sun et al. (2024), who detected wheat germination using NDVI time-series, the germination date for each field was derived as the x-axis position of μ - 2σ on the Gaussian fitted NDVI curve (see Figure 4a).

Flowering date retrieval

Zhou et al. (2023) showed that the peak of the GRVI time-series corresponds to the onset of wheat heading. Building on this, and given that flowering in wheat typically begins about 4 to 5 days after the completion of heading, which itself usually lasts 3–5 days (Large, 1954; Miller, 1992; Alabama Cooperative Extension System, 2022), the flowering date for each field was retrieved as the x-axis position of μ + 9 days on the Gaussian fitted GRVI curve (see Figure 4b), where μ represents the peak.

Senescence date retrieval

Anderegg et al. (2020) identified senescence onset on the increasing limb of the PSRI time-series as the point at which PSRI had increased by 20% of its amplitude from the trough, following the approach first applied to NDVI by Christopher et al. (2014) and subsequently adopted in many wheat senescence studies (Christopher et al., 2016; 2018; Kitonyo et al., 2017; Alhabbar et al., 2018). Building on this, and given that the fitted PSRI time-series in the present study was normalized between 0 and 1 (with a minimum of 0 and an amplitude of 1), the senescence date for each field was defined as the x-axis position corresponding to a normalized PSRI value of 0.20 on the Gaussian fitted PSRI curve (see Figure 4c), representing a 20% increase from the trough.

Harvesting date retrieval

Guided by the methodology of Pan et al. (2015), Sun et al. (2024), and El Imanni et al. (2022), who identified wheat harvesting using a NDVI time-series, the harvesting date for each field was retrieved as the x-axis position of μ + 2σ on the Gaussian fitted NDVI curve (see Figure 4a).

Cycle length estimation

The cycle length for each field was calculated as the difference between the NDVI-derived harvesting and germination dates, expressed as (μ + 2σ) - (μ - 2σ). Because cycle lengths were not extracted directly from the NDVI time-series but instead derived from NDVI-based retrievals, their agreement with the in situ observations provides an additional indicator of the reliability of the germination and harvesting date retrievals from Sentinel-2 data.

Retrieval of wheat phenology stages from Sentinel-1 data

In this section, we followed the same workflow described in Section 2.3 to retrieve wheat phenology stages using time-series of VV-VH derived from Sentinel-1 imagery. VV-VH reduces the influence of background effects, such as soil moisture and surface roughness, which tend to affect both VV and VH similarly (Ling et al., 2022). Therefore, this metric enhances the vegetation signal and emphasizes vegetation structure-related changes more effectively than either VV or VH alone, making it suitable for detecting canopy development and monitoring phenology.

As for the indices derived from Sentinel-2 imagery, VV-VH was generated for each Sentinel-1 image acquisition. For each image, the VV-VH values of all pixels within each delineated wheat field were averaged to construct field-level time-series profiles. Accordingly, a series of processing steps, closely aligned with those described in Section 2.3, were applied to the VV-VH time-series of each field in order to retrieve wheat phenology stages. First, each curve was normalized using the min-max method (Equation 4) (Vafaei et al., 2022). Then, the normalized data were smoothed using the Savitzky-Golay filtering technique (Equation 5) (Savitzky and Golay, 1964). Unlike the single Gaussian function applied in Section 2.3, the smoothed VV-VH time-series data were fitted using a two-component Gaussian function, expressed as follows (Kikawa et al., 2015): where f t is the fitted VV-VH value at time t; A 1 and A 2 are the amplitudes; μ 1 and μ 2 are the means; and σ 1 and σ 2 are the standard deviations of the first and second Gaussian components, respectively. The two-component Gaussian model was adopted for the following reasons: (1) visual inspection of the raw VV-VH time-series revealed a bimodal seasonal pattern, necessitating two components to capture the full temporal profile (see Figure 5); and (2) unlike the Sentinel-2 analysis in Section 2.3, which relied on three vegetation indices, the use of VV-VH data alone to retrieve all phenology stages required a more flexible model with additional curve parameters. Furthermore, the two-component Gaussian model was selected following a comparative evaluation of five candidate functions (single Gaussian, two-component Gaussian, polynomial, raised cosine, and Voigt) fitted to the VV-VH profiles. Across all wheat fields, the two-component Gaussian achieved the best fit quality, and was therefore adopted.

As in Section 2.3, curve-based rules were developed based on insights from previous wheat phenology studies (Nasrallah et al., 2019; Vavlas et al., 2020; Flores et al., 2025) to retrieve phenology stages from the fitted VV-VH curves. These studies employed fitted time-series of either VV-VH in decibels or VH/VV ratios in linear scale (using multiple Gaussian functions or logistic curves). Such SAR-based metrics are preferable for developing curve-based rules because they generally exhibit vegetation-driven seasonal trajectories and are less affected by soil moisture or irrigation events than single-polarization backscatter (VV or VH), thereby providing more consistent signals across environments (Vreugdenhil et al., 2020). Since the present study employed the VV-VH metric, it is worth noting that VV-VH expressed in decibels is mathematically equivalent to the logarithm of the VV/VH ratio in linear scale, and therefore represents the inverted form of the VH/VV ratio (i.e., VH/VV captures the same structural signal of VV-VH but with inverted polarity). This equivalence makes the transfer of established insights both reliable and feasible. The resulting rules were formulated as x-axis positions on the fitted VV-VH curves derived from the two-component Gaussian model parameters (A 1, A 2, μ 1, μ 2, σ 1, and σ 2; see Equation 7) and were applied individually to derive phenology stages for each field. While these rules originated from site-specific studies, their reliance on relative curve positions and the similarity in curve patterns support their adaptation to the current study context. The following subsections describe how each phenology stage was retrieved from Sentinel-1 data.

Germination date retrieval

Flores et al. (2025) identified wheat germination at the left-hand trough of the fitted VH/VV time-series. Building on this insight, and given the equivalence between VV-VH and VH/VV, the germination date for each field was defined as the x-axis position of μ 1 on the two-component Gaussian fitted VV-VH curve (see Figure 6), corresponding to the peak of the first Gaussian component.

Flowering date retrieval

Based on the work of Vavlas et al. (2020), who divided the VH/VV time-series into two phases, with the first fitted logistic curve (the left-hand portion of the seasonal VH/VV trajectory) representing the wheat growth phase used to identify flowering timing, the flowering date for each field was retrieved as the x-axis position of μ 1 + 2σ 1 on the two-component Gaussian fitted VV-VH curve (see Figure 6).

Senescence date retrieval

Vavlas et al. (2020) detected senescence onset by interpreting the second fitted logistic curve to the VH/VV time-series (the right-hand portion of the seasonal trajectory) as the wheat maturation phase and they associated the onset of this phase with the start of wheat senescence. Guided by this method, and given the equivalence between VV-VH and VH/VV, the senescence date for each field was derived as the x-axis position of μ 2 - 2σ 2 on the two-component Gaussian fitted VV-VH curve (see Figure 6).

Harvesting date retrieval

Nasrallah et al. (2019) detected wheat harvesting at the right-hand peak of the fitted VV-VH time-series. Given that the present study employed the VV-VH metric consistent with Nasrallah et al. (2019), the harvesting date for each field was retrieved as the x-axis position of μ 2 on the two-component Gaussian fitted VV-VH curve (see Figure 6), corresponding to the peak of the second Gaussian component.

Cycle length estimation

The cycle length for each field was calculated as the difference between the VV-VH-derived harvesting and germination dates, defined as (μ 2 - μ 1).

Accuracy assessment of retrieved wheat phenology stages

To evaluate the accuracy of wheat phenology stages retrieved independently from each of the Sentinel-2 and Sentinel-1 datasets, the derived DOY for each phenology stage was validated against the in situ data collected over the 2019–2020 and 2020–2021 growing seasons. For each phenology stage, one scatter plot for the Sentinel-2-based retrievals and another for the Sentinel-1-based retrievals were generated, each combining data from both growing seasons into a single regression analysis. This combined-season approach was not applied to flowering and senescence stages, which were retrieved for a single season only because the corresponding in situ observations were available exclusively for the 2020–2021 season, as noted in Section 2.2.3. The combined-season analysis increased sample size and improved the robustness of the statistical evaluation, and was justified by the similar variability patterns in crop growth observed in the in situ data across the two seasons (see Table 1). The accuracy assessment employed the coefficient of determination (R 2) (Equation 8), root mean square error (RMSE) (Equation 9), relative RMSE (rRMSE) (Equation 10), and bias (Equation 11), as defined below: where n is the number of observations; x i and y i are the observed and derived values, respectively; and x ̄ and ȳ are their corresponding means. Additional details regarding these statistical metrics can be found in (Krause et al., 2005; Ma et al., 2011; Maupin et al., 2019; Mabasa et al., 2020). For DOY values that straddle the end of the calendar year for the same phenology stage, such as those related to germination stage, a continuous DOY counter spanning the two calendar years of each growing season was used for the computations of the accuracy metrics defined above, such that dates occurring after 31 December continued to increment past DOY 365 rather than resetting. This ensured that date differences, and the accuracy metrics, were computed correctly.

The accuracy assessment aimed both to explore the performance of each satellite dataset and to examine the applicability of the curve-based rules under dryland conditions. Examining the applicability of the curve-based rules is particularly critical because, although these rules were developed based on established insights from previous studies, with the specific study underpinning each rule identified in Sections 2.3 and 2.4, dryland environments pose unique climatic and management challenges that may influence their performance or alter the temporal patterns of satellite data. A curve-based rule was considered inapplicable when its retrieval error, specifically rRMSE, reached or exceeded 15%. In such cases, additional fine-tuning was undertaken to improve the predictive performance of the rule and reduce the associated error.

Alignment of the curve-based rules with the in situ mean phenology dates

The min-max normalized NDVI, GRVI, PSRI, and VV-VH time-series, averaged across all active wheat fields during the corresponding growing season, together with their corresponding Gaussian-fitted curves, are illustrated in Figure 7. For each index, the dates of relevant phenology stages observed in situ (computed as the average across all fields and referred to as the “mean date”), along with the rule-based dates derived from the fitted mean curves, were plotted to visualize their alignment and to initially assess the ability of the rules to capture the mean dates observed in situ. To avoid visual clutter, we displayed only means (excluding standard deviations and 95% confidence bands).

The germination mean dates observed in situ were well identified by the rules of both NDVI and VV-VH in both seasons. Applying the NDVI rule for detection of germination resulted in delays of 1 week in 2019–2020 and between one and 2 weeks in 2020–2021 compared to the in situ observations (see Figures 7a,b). In contrast, the VV-VH rule identified germination between one and 2 weeks prior to the in situ information in both seasons (see Figures 7e,f). For flowering detection, the GRVI rule showed superior performance compared to the VV-VH rule. The GRVI rule predicted flowering prior to the flowering mean date observed in situ between one and 2 weeks in 2020–2021 (see Figure 7c). The VV-VH rule placed flowering between two and 3 weeks earlier than the in situ observation (see Figure 7f). The senescence mean date observed in situ was well detected by the rules of both PSRI and VV-VH. The PSRI rule predicted senescence earlier than the in situ observation between one and 2 weeks in 2020–2021 (see Figure 7d). Applying the VV-VH rule yielded earlier detection of senescence relative to the in situ mean date by 1 week in 2020–2021 (see Figure 7f). The harvesting stage proved more difficult to detect, with the NDVI rule placing harvesting later than the harvesting mean dates observed in situ by three to 4 weeks in 2019–2020 and 2 weeks in 2020–2021 (see Figures 7a,b). The VV-VH rule performed worst, with harvesting time lagging behind the in situ observations by more than a month in 2019–2020 and almost a month in 2020–2021 (see Figures 7e,f).

Earlier detection of a phenology stage date compared to the mean date observed in situ suggested that the corresponding rule responded not to the timing of the phenology stage itself, but rather to a precursor condition that preceded and was likely associated with that stage. In contrast, delayed detection indicated that the rule responded directly to the phenology stage timing. Although the analysis in this section focused on averages across all fields, these visual comparisons provided an informative narrative interpretation and an initial indication that the curve-based rules are broadly applicable in the dryland agricultural setting examined in this study: though clear stage- and index-specific differences were evident. The subsequent Sections 3.2-3.6 present the full quantitative evaluation, where field-level scatterplots and error metrics are used to confirm which rule-index pairings are most suitable and where refinements may be required.

Germination date retrieval from Sentinel-2 and from Sentinel-1 data

Germination dates were retrieved from NDVI and VV-VH time-series using the respective curve-based rules across the 68 and 73 active wheat fields during the 2019–2020 and 2020–2021 growing seasons, respectively. Figure 8 compares the in situ germination dates with those derived from the fitted NDVI (Figure 8a) and VV-VH (Figure 8b) time-series for individual wheat fields. For the NDVI-based germination dates (see Figure 8a), the regression line was shifted upward the 1:1 line, reflecting an image-predicted delay in germination (also see Figures 7a,b), with an RMSE of 9.03 days, an rRMSE of 2.64%, a bias of 7.34 days, and an R 2 of 0.83. The high R 2 value showed that the NDVI-based germination dates explained 83% of the field-to-field variability in the in situ germination dates. Furthermore, the positive bias suggests that the NDVI rule responded primarily to the signal of wheat germination, with the 7-day delay, reflecting the image availability and 5-day temporal sampling around the germination window. This delay likely reflects the nature of the optical measurement. NDVI responds to the contrast between near-infrared and red reflectance produced by photosynthetically active material, which becomes detectable only once emerging seedlings have developed sufficient leaf area to dominate the reflectance of the underlying soil. Germination is therefore registered by NDVI once the canopy has developed enough to alter the field-level optical signal.

The VV-VH-based germination dates (see Figure 8b) exhibited the opposite tendency to those derived from NDVI, with the regression line shifted below the 1:1 line, indicating early detection (also see Figures 7e,f), with an RMSE of 15.75 days, an rRMSE of 4.60%, a bias of −12.93 days, and a lower R 2 of 0.68. The R 2 value showed that the VV-VH retrievals still explained a substantial portion (nearly 70%) of the field-to-field variability in the in situ germination dates. However, the negative bias suggests that the VV-VH rule responded to surface or structural conditions occurring an average of 13 days before the in situ observed germination dates, rather than to the germination signal itself. Given this 13-day lead, which closely matches the sowing-to-germination interval for wheat (Wang et al., 2009) with an average of 12.2 days across 244 sites, the precursor conditions detected by the VV-VH rule could plausibly be linked to field-preparation processes occurring right before or shortly after sowing. Within this context, the VV-VH rule may respond to one of several possible factors: for instance, the sowing activity itself, soil disturbance, changes in surface roughness, moisture redistribution following sowing, or other early surface alterations. This behavior has also been documented in (Song and Wang, 2019; Nasirzadehdizaji et al., 2021), in which Sentinel-1 backscatter responds strongly to sowing operations and the accompanying field-preparation disturbances. We emphasize that this interpretation is presented as a hypothesis consistent with the available evidence discussed above, rather than as a confirmed mechanism, since concurrent in situ soil moisture and surface roughness measurements were not collected in this study.

Accuracy metrics computed separately for each growing season are provided in Supplementary Table S1, and confirm that agreement with the in situ observations is maintained within each individual season. The NDVI rule offered direct retrieval of the germination timing, whereas the VV-VH rule captured pre-germination field activity (might be field-preparation processes occurring around wheat sowing). It is essential to recognize this distinction when interpreting the results or applying the rules operationally. Overall, these findings highlight the complementary value of the NDVI and VV-VH rules for the operational retrieval of wheat germination dates and, more broadly, for assessing early-season phenology (including both germination and pre-germination) under dryland conditions.

Flowering date retrieval from Sentinel-2 and from Sentinel-1 data

Flowering dates were derived from GRVI and VV-VH time-series for 73 active wheat fields during the 2020–2021 growing season, and validated against the in situ flowering dates. The GRVI-based flowering dates (see Figure 9a) showed early retrieval relative to in situ observations (also see Figure 7c), yielding an RMSE of 9.63 days, an rRMSE of 10.46%, a bias of −9.28 days, and an R 2 of 0.55. The early detection of flowering dates obtained with the GRVI rule may be attributed to the design of the GRVI rule itself. It is the only rule in our framework that applies a fixed temporal shift (rather than a proportional or curve-derived percentage adjustment), specifically adding 9 days to the mean of the fitted GRVI curve (the GRVI-derived peak), as described in Section 2.3. This 9-day assumption originates from the notion that the wheat heading stage (corresponding to the GRVI peak as reported in (Zhou et al., 2023)) may begin and conclude within roughly 3–5 days, followed by the onset of flowering around 4 days later. However, this fixed 9-day shift does not appear to be suitable under the environmental conditions of the study site, which is evident from the observed negative bias, causing the rule to predict flowering earlier than the in situ dates. The results suggest that a larger shift would be more appropriate and would reduce the magnitude of the resulting negative bias. Despite that, the R 2 value obtained here indicated that 55% of the field-to-field variability in the in situ flowering dates was explained by the GRVI retrievals. Beyond the design of the GRVI rule, this behavior may reflect what the index itself measures. GRVI is computed from the contrast between green and red reflectance and therefore tracks the change in canopy pigment composition that accompanies the transition from vegetative growth into reproductive development. Since this pigment transition culminates at heading rather than at flowering, the resulting curve landmark precedes flowering, and the retrieval is displaced early by however much the assumed interval between the two falls short.

The VV-VH-based flowering dates (see Figure 9b) demonstrated the same overall tendency as the GRVI retrievals but with a more pronounced downward shift of the regression line from the 1:1 line. The associated accuracy metrics yielded an RMSE of 18.38 days, an rRMSE of 19.95%, a bias of −18.02 days, and an R 2 of 0.60. The early retrieval observed for the VV-VH-based flowering dates may be explained by the likelihood that the VV-VH rule responded to an earlier phase in wheat reproductive development, either booting or heading, rather than to flowering itself. Because the observed in situ dates for booting or heading were not available, it is difficult to determine exactly what the VV-VH rule responded to. However, Zhou et al. (2023) reported that the GRVI peak consistently corresponds to the wheat heading stage. Hence, we visually compared the VV-VH flowering rule in Figure 7f with the GRVI peak in Figure 7c and found that the VV-VH rule does not align with the GRVI peak. Instead, it occurs earlier, indicating that the VV-VH rule is not responding to heading but rather to an even earlier phenology phase. The most plausible stage is booting, which directly precedes heading in wheat development. During the booting stage, several physiological and structural dynamics occur that VV-VH backscatter may respond to, including rapid structural reorganization of the upper canopy, early spike swelling within the flag leaf sheath, and moisture-driven dielectric changes. This indicates that the current formulation of the VV-VH rule requires fine-tuning before it can be reliably used for wheat flowering retrieval. Such refinement is both expected and appropriate, given that the existing rule was adapted from (Vavlas et al., 2020) and was not specifically optimized for irrigated hot dryland conditions. Under dryland conditions, the GRVI rule provides the more dependable pathway for the operational retrieval of wheat flowering dates, while the VV-VH rule may still hold value but only after targeted refinement.

Senescence date retrieval from Sentinel-2 and from Sentinel-1 data

The dates of senescence were retrieved using the respective curve-based rules applied to the PSRI and VV-VH time-series for 73 active wheat fields during the 2020–2021 growing season. Figure 10 compares the in situ senescence dates with those derived from the fitted PSRI (Figure 10a) and VV-VH (Figure 10b) time-series for individual wheat fields. For the PSRI-based senescence dates (see Figure 10a), the regression line occurred below the 1:1 line, indicating early retrieval (also see Figure 7d), with an RMSE of 12.22 days, an rRMSE of 10.35%, a bias of −11.97 days, and an R 2 of 0.59. Likewise, the VV-VH-based senescence dates (see Figure 10b) followed a similar pattern but with a smaller downward shift of the regression line from the 1:1 line (also see Figure 7f), with an RMSE of 7.74 days, rRMSE of 6.56%, a bias of −7.11 days, and a higher R 2 of 0.63. Retrievals based on both PSRI and VV-VH were able to explain approximately comparable proportions of the field-to-field variability in the in situ senescence dates.

The early retrieval obtained by both rules can be explained by the criteria used by farmers to assign senescence dates in situ and the specific senescence timing detected by the applied image-based rules. Unlike other phenology stages considered in this study, senescence is inherently a gradual and progressive process, characterized by a slow canopy dry-down in which chlorophyll is progressively lost, nitrogen and carbohydrates are remobilized into the grain, and leaf color transitions from green to yellow and ultimately brown over an extended period. This process spans several weeks and, depending on environmental conditions, may extend over longer periods (Borrill et al., 2019), making senescence less temporally discrete than other stages that start and conclude more rapidly. In this context, it is important to distinguish between the senescence timing that the rules were designed to detect and the field-based criterion used to assign senescence dates. In fact, the PSRI and VV-VH rules applied here were explicitly designed to capture the onset of senescence, as reported in the studies from which these rules were adopted (Anderegg et al., 2020; Vavlas et al., 2020), and as further supported by the physical interpretation of PSRI and VV-VH at the corresponding rule locations along the time series. Specifically, the PSRI rule location reflects the initiation of chlorophyll degradation and the emerging dominance of carotenoids that together mark the senescence onset before widespread yellowing. In parallel, the location of the VV-VH rule reflects an increase in VV relative to VH, indicating a reduction in depolarization associated with dense green vegetation and a transition toward simpler and drier canopy structures characteristic of early senescence. In contrast, farm records confirmed that in situ senescence dates were assigned when a substantial portion of the wheat canopy exhibited yellowing, corresponding to the midpoint rather than the onset of senescence. Consequently, the observed early detection reflects a systematic difference between the senescence onset captured by the PSRI and VV-VH rules and the field-based criterion used for in situ senescence observations, rather than a misidentification of unrelated pigment or structural precursor changes. This interpretation is supported by two related points. First, the magnitude of the associated biases (7 and 11 days for the VV-VH and PSRI rules, respectively) is consistent with the expected temporal separation between senescence onset and an intermediate stage within a gradual process spanning several weeks. Second, the small difference in bias between the PSRI and VV-VH rules (around 4 days) indicates that both rules responded to the same event, with the 4-day difference being attributed to variations in temporal sampling and signal response between Sentinel-2-derived PSRI and Sentinel-1-derived VV-VH, rather than to the detection of different events.

Importantly, both rules were able to respond to the wheat senescence stage itself (specifically its onset) rather than to a prior phenology stage or unrelated pigment or structural precursor changes. Accordingly, the retrieved dates should be interpreted as representing the onset of wheat senescence specifically rather than later senescence timing. Overall, these results highlight the applicability of both PSRI and VV-VH rules for the operational retrieval of wheat senescence onset dates under dryland conditions.

Harvesting date retrieval from Sentinel-2 and from Sentinel-1 data

Harvesting dates were derived for 141 active wheat fields during the two growing seasons from NDVI and VV-VH time-series. Figure 11 shows the relationships between the in situ harvesting dates and those derived from the fitted NDVI (Figure 11a) and VV-VH (Figure 11b) time-series for individual wheat fields. The NDVI-based harvesting dates (see Figure 11a) revealed a delayed detection (Figures 7a,b), with the regression occurring above the 1:1 line. The associated accuracy metrics produced an RMSE of 21.26 days, an rRMSE of 14.17%, a bias of 20.60 days, and an R 2 of 0.34. The VV-VH-based harvesting dates (see Figure 11b) exhibited an even larger detection delay (Figures 7e,f) reflected in the accuracy metrics, including an RMSE of 34.75 days, an rRMSE of 23.70%, a bias of 34.36 days, and a lower R 2 of 0.24. Both rules demonstrated low R 2 values, indicating that the retrieved harvesting dates explained only a limited proportion of the field-to-field variability in the in situ harvesting dates.

While smaller retrieval errors and biases were typically expected under the assumption that pre- and post-harvest NDVI and VV-VH signal changes would be sharp and detectable within the revisit frequency of Sentinel observations, the derived harvesting dates from both rules exhibited larger deviations and systematic delays relative to the in situ dates. These delays can be attributed to several factors. First, the NDVI and VV-VH rules applied here were designed to detect the point at which wheat harvesting is completed across the entire field, whereas the recorded in situ harvesting date reflects the start rather than the completion of the harvesting process. Given these design considerations, and given that harvesting is rarely instantaneous at the field scale and often extends over several days, harvesting dates derived from both rules are thus expected to be delayed relative to in situ dates. Second, NDVI and VV-VH are not directly sensitive to the harvesting operation itself but instead respond to physiological and structural changes that occur progressively during senescence and maturity. This results in smooth temporal transitions rather than abrupt changes at harvest, meaning that harvesting date detection is not necessarily constrained by the revisit frequency of Sentinel observations and can therefore exhibit larger deviations. Third, the use of field-averaged NDVI and VV-VH values further smooths spatial variability during harvesting, as residual vegetation can maintain non-minimal index values even when most of the field has been harvested, thereby delaying the point at which the field is identified as fully harvested, which is the condition targeted by both rules. Although these factors affect both rules, their combined influence may be expressed differently in NDVI and VV-VH, leading to larger delays in the VV-VH retrievals compared to the NDVI retrievals (20.60 and 34.36 days for the NDVI and VV-VH rules, respectively), potentially due to differences in the sensitivity of VV-VH to harvesting conditions relative to NDVI. In addition to this, gaps in the effective acquisition sequence of Sentinel-1 may also exacerbate retrieval delays. Referring to Table 1, in situ harvesting during the 2019–2020 season occurred between DOY 142 and 152 (21–31 May), during which Sentinel-1 exhibited a 12-day gap in late May, while Sentinel-2 acquisitions remained at 5-day intervals (see Figures 2, 3). Similarly, during the 2020–2021 season, in situ harvesting occurred between DOY 136 and 155 (16 May-4 June), where Sentinel-1 again showed a 12-day gap in mid-May, whereas Sentinel-2 coverage remained temporally dense (see Figures 2, 3). Such gaps can further delay harvesting date detection by the VV-VH rule by limiting observations immediately following in situ harvest completion. The difference in delay between the sensors reflects what each measures as the crop is removed. By harvest, the canopy has lost most of its chlorophyll through senescence, so NDVI trajectory has fallen close to its seasonal minimum before cutting. As NDVI is driven by red and near-infrared contrast associated with photosynthetically active material, removal of already-senesced crop produces only a small change relative to the seasonal amplitude. Volume backscattering behaves differently, as radar responds to the geometry and water content of the material present rather than its pigment. Standing stubble and crop residue continue to contribute to volume scattering after cutting, since dry stems retain the structural arrangement to which the signal responds, and the VV-VH retrieval may be thus displaced further from the recorded harvesting date.

Accuracy metrics computed separately for each growing season are provided in Supplementary Table S1, and confirm that agreement with the in situ observations is maintained within each individual season. The VV-VH rule showed limited reliability in its current form, as indicated by a larger bias and an rRMSE exceeding the 15% threshold at approximately 24%, highlighting the need for further rule refinement to improve harvesting date retrieval. Such refinement is expected and appropriate, given that the VV-VH harvesting rule was adapted from (Nasrallah et al., 2019) and was not specifically optimized for the local conditions considered here, similar to the VV-VH flowering rule discussed in Section 3.3. It is also noteworthy that, despite the observed delays, both rules here directly targeted the harvesting stage itself rather than earlier phenology stage or indirect spectral or structural precursors, with the delayed detection likely attributable to the factors discussed above. Overall, these findings support the suitability of the NDVI rule for the operational retrieval of wheat harvesting dates under dryland conditions, while suggesting that the VV-VH rule would benefit from targeted fine-tuning to reduce the systematic offset observed between the rule-derived and in situ dates.

Cycle length retrieval from Sentinel-2 and from Sentinel-1 data

Cycle lengths were calculated as the difference between the individually retrieved harvesting and germination dates for 141 active wheat fields across the two growing seasons, and validated against the in situ cycle lengths. For the NDVI-based cycle lengths (see Figure 12a), the regression line was shifted upward the 1:1 line, showing overestimation compared to the in situ observations with an RMSE of 16.32 days, an rRMSE of 9.69%, a bias of 13.15 days, and an R 2 of 0.65. The VV-VH-based cycle lengths (see Figure 12b) exhibited an even larger magnitude of overestimation, with an RMSE of 44.38 days, an rRMSE of 26.22%, a bias of 43.73 days, and a similar R 2 of 0.65. Despite the differences in absolute accuracy, cycle lengths based on both NDVI and VV-VH explained an equivalent proportion (65%) of the field-to-field variability in the in situ cycle lengths. The reported errors and biases in cycle length estimates derived from both NDVI and VV-VH were influenced by the accuracy of the germination and harvesting dates derived from NDVI and VV-VH, given that cycle length was calculated as the difference between those two phenology stages. Accuracy metrics computed separately for each growing season are provided in Supplementary Table S1, and confirm that agreement with the in situ observations is maintained within each individual season. Taken together, these results highlight the suitability of the NDVI-based approach for the estimation of wheat cycle lengths under dryland conditions, while also indicating that the VV-VH approach could become a viable option if the underlying rules (particularly the harvesting rule) are further refined and optimized for local conditions.

Refinement of curve-based rules with low retrieval accuracy

Rules were applied herein to time-series of optical vegetation indices and volume backscattering produced retrievals with varying levels of performance across wheat phenology stages. Retrievals from optical indices showed rRMSE ranging from 2.64% to 14.17%, while those from volume backscattering ranged from 4.60% to 26.22%. Since all these rules were developed based on insights from prior studies without any adjustment specific to our local conditions, there is clear scope for refining rules that showed high retrieval error (those with rRMSE exceeding the 15% threshold, as defined in Section 2.5), which are the VV-VH flowering and harvesting rules. Refining these two rules involves adjusting their current form to reduce the offsets observed in Figures 7e,f between the rule lines and the in situ references. Therefore, two adjustments are proposed. For flowering, the current VV-VH rule μ 1 + 2σ 1 (see Figure 6; adapted from (Vavlas et al., 2020)) predicted flowering between two and 3 weeks prior to the flowering mean date observed in situ (see Section 3.1; see Figure 7f). Shifting this rule rightward to the x-axis position of μ 1 + 3σ 1 is expected to reduce the current offset between the rule line and the in situ reference. For harvesting, the current VV-VH rule μ 2 (see Figure 6; adapted from (Nasrallah et al., 2019)) placed harvesting later than the in situ observations by more than a month in 2019–2020 and almost a month in 2020–2021 (see Section 3.1; see Figures 7e,f). Shifting this rule leftward to μ 2 - σ 2 should correct the current offset. Because cycle length estimation depends on both germination and harvesting dates, refining the VV-VH harvesting rule should likewise improve the accuracy of VV-VH-based cycle length estimates. Figure 13 illustrates the fine-tuned flowering and harvesting rules, which now place both stages less than a week later than the in situ observations, thereby substantially reducing the original offsets shown in Figures 7e,f. Retrievals based on these fine-tuned rules were obtained and then quantitatively evaluated against the corresponding in situ observations for all wheat fields (see Figure 14), yielding lower errors than the original rules (see Figures 9b, 11b). Furthermore, the positive biases associated with retrievals based on the fine-tuned rules (see Figures 14a,b) suggest that both rules responded primarily to the timing of the phenology stages they were refined to retrieve, rather than to precursor conditions preceding those stages. It is noteworthy that this evaluation used the same in situ observations from which the offsets were identified, and therefore demonstrates that the flowering and harvesting errors are systematic and geometrically correctable offsets rather than random error, but does not by itself constitute an independent validation of the refined rules. Two aspects of the procedure should nonetheless be noted. First, each rule was adjusted by a single σ of the fitted curve rather than by a numerically optimized value, so the refinement is a principled geometric correction of a systematic offset rather than a fit of free parameters to the observations. Second, the offsets were identified from the mean curve averaged across all fields and the in situ mean date (see Figure 13), whereas the resulting rules were applied and evaluated at the level of each individual field.

Other rules that produced retrievals with rRMSE below the 15% threshold were not adjusted, as an rRMSE below this threshold is reasonably considered indicative of good agreement in a satellite-based validation context. That rationale is based on the inherent uncertainties in both remote sensing and in situ observations, as well as within-farm variability in management practices and irrigation schedules, which may introduce spectral noise into satellite signals and reduce retrieval accuracy (Bégué et al., 2018). In combination, these factors highlight the need to interpret phenology retrieval results in the context of satellite limitations and real world-conditions.

With the performed refinements to the VV-VH flowering and harvesting rules, and given the already solid performance of the VV-VH germination and senescence rules, the VV-VH data has the potential to become a viable option for wheat phenology detection alongside the NDVI, GRVI, and PSRI data, subject to independent validation of the refined rules on data not used to derive them. This is crucial for operational wheat phenology monitoring in cloud-prone regions, where Sentinel-1 data remains a vital asset due to its all-weather imaging capabilities (Frison et al., 2018).

Sensing characteristics of Sentinel-2 and Sentinel-1 and their influence on retrieval accuracy

The retrieval accuracies obtained across the five phenology stages reflect a consistent pattern that follows from the different nature of the two observing systems. Sentinel-2 is a passive optical instrument measuring reflected solar radiation (Spoto et al., 2012), and the vegetation indices derived from it respond to canopy pigment composition and greenness. Sentinel-1 is an active radar system (Attema et al., 2009), and the volume backscattering derived from it responds to the geometry, arrangement, and water content of the canopy (Vreugdenhil et al., 2020; Ling et al., 2022). Each phenology stage examined here is expressed through changes in one or both of these properties, and the accuracy achieved by a given sensor at a given stage depends on how directly the transition being retrieved is expressed in the property that sensor measures.

Such correspondence accounts for the stage-dependent behavior observed in the Results. Germination and flowering are expressed most immediately through the development and pigmentation of photosynthetic material, and the optical retrievals track them more closely, whereas the radar retrievals respond to associated physical changes that precede them, namely, the surface disturbance accompanying sowing and the structural reorganization of the canopy during booting. Senescence, by contrast, is a process in which pigment change and structural change proceed together, as chlorophyll degradation is accompanied by canopy drying and the loss of the dense, randomly oriented green structures that produce depolarization (Vreugdenhil et al., 2020). Both sensors therefore respond to genuine senescence dynamics, and both retrievals show comparable agreement with the in situ observations. Harvesting is expressed least directly in either property, since neither pigment composition nor canopy structure changes abruptly at the moment of cutting, and the residual structural signal from stubble and crop residue persists after the optical signal has already reached its seasonal minimum.

Implications, limitations, and potential future work

The findings outlined in this study have valuable implications for operational wheat phenology detection. By applying simple curve-based rules to freely available Sentinel-2 and Sentinel-1 data, each analyzed independently, we showed that multiple phenology stages can be retrieved with reasonable accuracy in a dryland farming region, where satellite time-series interpretation is typically complicated by strong soil-vegetation signal mixing under such environments (Smith et al., 2019). This highlights the potential of a practical and scalable method, which can be applied on wheat fields in any season once normalized, smoothed, and Gaussian-fitted field-level time-series data obtained from satellite imagery for the wheat fields are available in advance. The simple approach used herein substitutes for labor-intensive ground-based phenology observations that are often sparse and infeasible at scale. It enables reliable and low data-requirement extraction of phenology information essential as input for crop simulation models (Meng et al., 2009). Such capabilities offer a promising future work direction that could facilitate crop model calibration under limited ground-based observations and support large-scale crop-water-atmosphere modeling. Beyond modeling, phenology data also improves crop type classification and condition assessment (Chu et al., 2016). Another strength of the study lies in the inclusion of two different data sources (Sentinel-2 and Sentinel-1), revealing their respective capabilities across five phenology stages. Moreover, the analysis was grounded in extensive field observations from two growing seasons, ensuring the evaluation reflected real production conditions rather than controlled experiments, thereby underscoring the robustness and relevance of the findings.

Several limitations and potential future works should be acknowledged. The analysis was limited to a single large farm, meaning that the generalizability and transferability of the rules to other dryland systems require careful consideration. Future research should extend the analysis and validation across diverse dryland environments and management practices to further test the broad applicability of the rules under more conditions that may exist in drylands. Establishing this broader applicability will also require a systematic sensitivity analysis characterizing the responsiveness of the retrieved dates to the rule parameters, which was beyond the scope of the present single-site evaluation. In addition, the present evaluation was based on two growing seasons and a single wheat cultivar. While the two seasons differed appreciably in sowing schedule and growth duration, and while the use of a single cultivar isolated management- and environment-driven phenology variability from varietal effects, evaluation across additional cultivars and growing seasons remains necessary to establish the transferability of the rules to phenology expressions and seasonal conditions beyond those examined here. Furthermore, the evaluation conducted herein was retrospective, with rules applied after the growing season, restricting the current approach to post-season analysis. Future efforts should prioritize transitioning the current approach from retrospective to near- or real-time, which would enable operational in-season monitoring. Within this context, the utility of the current approach would be enhanced and the retrieved phenology information would support in-season agronomic decision-making, such as irrigation scheduling, fertilizer application, and pest control (Baghdadi et al., 2009; Fieuzal et al., 2013; Mandal et al., 2018). Moreover, the refined VV-VH flowering and harvesting rules were evaluated on the same in situ observations used to identify the offsets they correct. Independent validation of these refined rules on data not used to derive the adjustments is therefore required before they are applied operationally, and represents a necessary direction for future work. Finally, the satellite datasets themselves impose limitations. For Sentinel-1, volume backscattering initially revealed limited sensitivity for flowering and harvesting stages, likely requiring a preliminary visual alignment with in situ observations as demonstrated in Figure 7. In addition, in the center-pivot irrigated setting studied here, irrigation events can transiently alter surface and canopy moisture and thereby introduce short-term variations in the backscatter signal, and to a lesser extent in the optical signal, that are unrelated to phenology development. Such effects are most likely to influence retrieval during the low-biomass stages, when the canopy contribution to the signal is small, and their explicit characterization represents a useful direction for improving rule-based retrieval in irrigated dryland systems. For Sentinel-2, the optical indices depend on cloud-free observations, which may constrain accuracy in cloudy regions (Misra et al., 2020), unless supported by gap-filling strategies (Kandasamy et al., 2013), or through integration with additional sensors (Claverie et al., 2018). A further consideration for Sentinel-2 is that the approach requires three vegetation indices rather than a single one. We examined whether it could be simplified by retrieving all stages from NDVI alone, applying structural analogues of the GRVI and PSRI rules to the fitted NDVI curves. This produced substantially larger errors and systematic offsets for flowering and senescence (see Supplementary Table S2; Supplementary Figures S1, S2), indicating that the reliance on stage-specific indices is necessary rather than incidental.

The employed retrieval approach is based on the geometry of Gaussian curves fitted to normalized and smoothed time series rather than absolute signal values. Furthermore, the fitted parameters used to design the rules are principally determined by the high-amplitude portion of the seasonal trajectory, when canopy cover is near-complete and soil contribution to the observed signal is minimal. Consequently, the approach is relatively insensitive to constant shifts caused by soil background, although this effect differs between the two sensors. For Sentinel-1, the use of the VV-VH backscatter difference further minimizes soil background influence. For Sentinel-2, however, soil background influence is not fully eliminated, specifically during low-biomass periods, which coincide with the germination and harvesting rule positions in the tails of the fitted curves. Thus, residual soil effects may contribute to the retrieval errors reported for these stages. Future studies should evaluate soil-adjusted vegetation indices (e.g., SAVI and OSAVI) against the indices employed here to determine whether they further improve phenology retrieval in bright sandy loam soils under dryland conditions.

This study evaluated the feasibility of using freely available Sentinel-2 optical vegetation indices (NDVI, GRVI, and PSRI) and Sentinel-1 volume backscattering, combined with curve-based rules, to retrieve wheat phenology stages for a large number of fields cultivated in a dryland setting. Five stages (germination, flowering, senescence, harvesting, and cycle length) were retrieved from each satellite separately and validated against extensive in situ observations from two growing seasons. Results showed that rules applied to NDVI, GRVI, and PSRI produced reasonably accurate retrievals across all phenology stages without requiring adjustment. Rules applied to VV-VH showed improved performance across all phenology stages following targeted refinement of the flowering and harvesting rules, which initially showed high retrieval errors, unlike the VV-VH germination and senescence rules that performed reasonably well without modification. This improvement was assessed on the same in situ data used to identify the refinements and therefore requires independent validation. Importantly, the employed rules responded directly to the phenology stages they were designed to retrieve. A notable exception was the VV-VH germination rule, which responded to field-preparation processes around sowing, indicating complementary value alongside the NDVI germination rule for broadly assessing early-season phenology. A critical insight is that discrepancies between satellite retrievals and field observations sometimes arise from targeting different moments within the same phenology stage. Specifically, rules may detect the onset or completion of a stage, while field observations often represent an intermediate, more visually apparent moment in the same stage progression. This was particularly notable for the senescence and harvesting stages. Overall, the approach utilized herein offers a practical, scalable, and resource-efficient solution for retrieving wheat phenology stages in the data-scarce and operationally constrained dryland farming region examined in this study, reducing dependence on labor-intensive and spatially limited ground-based observations. The resulting phenology information can support crop modeling, improve crop classification, and inform in-season agricultural management decisions (while recognizing that such use would require transitioning the approach from retrospective to near- or real-time application). It is also crucial to note that the approach was demonstrated under a specific set of conditions, i.e., sandy-loam soils, center-pivot irrigation, and a hot desert (BWh) climate. Its applicability to other soil textures, irrigation systems, and climatic regimes should therefore be established through targeted evaluation before the rules are applied operationally in those settings. Collectively, this work establishes a foundation for advancing operational wheat phenology monitoring in dryland agriculture, underscoring the strengths of applying simple curve-based rules to Sentinel-2 and Sentinel-1 data, each treated independently, and highlighting opportunities for expansion and integration into broader agricultural decision-support systems when the approach is tested over diverse dryland environments.