Prior Season Crop Type Masks for Winter Wheat Yield Forecasting: A US Case Study

Abstract

Monitoring and forecasting crop yields is a critical component of understanding and better addressing global food security challenges. Detailed spatial information on crop-type distribution is fundamental for in-season crop condition monitoring and yields forecasting over large agricultural areas, as it enables the extraction of crop-specific signals. Yet, the availability of such data within the growing season is often limited. Within this context, this study seeks to develop a practical approach to extract a crop-specific signal for yield forecasting in cases where crop rotations are prevalent, and detailed in-season information on crop type distribution is not available. We investigated the possibility of accurately forecasting winter wheat yields by using a counter-intuitive approach, which coarsens the spatial resolution of out-of-date detailed winter wheat masks and uses them in combination with easily accessibly coarse spatial resolution remotely sensed time series data. The main idea is to explore an optimal spatial resolution at which crop type changes will be negligible due to crop rotation (so a previous seasons’ mask, which is more readily available can be used) and an informative signal can be extracted, so it can be correlated to crop yields. The study was carried out in the United States of America (USA) and utilized multiple years of NASA Moderate Resolution Imaging Spectroradiometer (MODIS) data, US Department of Agriculture (USDA) National Agricultural Statistics Service (NASS) detailed wheat masks, and a regression-based winter wheat yield model. The results indicate that, in places where crop rotations were prevalent, coarsening the spatial scale of a crop type mask from the previous season resulted in a constant per-pixel wheat proportion over multiple seasons. This enables the consistent extraction of a crop-specific vegetation index time series that can be used for in-season monitoring and yield estimation over multiple years using a single mask. In the case of the USA, using a moderate resolution crop type mask from a previous season aggregated to 5 km resolution, resulted in a 0.7% tradeoff in accuracy relative to the control case where annually-updated detailed crop-type masks were available. These findings suggest that when detailed in-season data is not available, winter wheat yield can be accurately forecasted (within 10%) prior to harvest using a single, prior season crop mask and coarse resolution Normalized Difference Vegetation Index (NDVI) time series data.

1. Introduction

In recent years, there has been a dramatic increase in the demand for up-to-date, reliable global agricultural information [1,2,3,4,5]. The issue of food security has rapidly risen to the top of government agendas around the world following the food price spikes in August 2007 and November 2010. Consequently, there is a growing recognition that more timely and accurate information on global crop production is indispensable for combating the growing stress on the world’s crop production, stabilizing food prices, developing effective agricultural policies, and coordinating responses to regional food shortages [6].

Remotely sensed data from space offers a practical means for generating such information, as they provide consistent, global, timely, cost-effective, and synoptic information on crop condition and distribution. Their utility for crop yield forecasting has long been recognized and demonstrated across a wide range of scales and geographic regions [7,8,9,10,11] Nevertheless, it is widely acknowledged that remotely sensed data could be better utilized within current operational monitoring systems, so there is a critical need for research focused on developing practical robust methods for agricultural monitoring and specifically for crop production forecasting using available and accessible data [1,12,13].

Stratifying a region into different crop types—a process commonly termed as crop masking—is an important step in developing Earth Observations (EO)-based yield forecasting models [14,15]. Such masks enable the isolation of the remotely sensed (RS) crop-specific signal throughout the growing season, reducing the noise on the signal from other land cover or crop types. A variety of crop masking methods used for isolating overall cropland signals as well as crop-specific signals for the purpose of yield forecasting. Atzberger [16] provides an extensive review of existing remote-sensing-based agriculture monitoring systems. However, one of the difficulties in monitoring and forecasting crop yields using RS imagery is the availability of timely seasonal detailed crop type masks that can be used to identify the crop of interest prior to the end of the growing season. Often, a general cropland mask is used to distinguish cropped areas from other land use types, rather than a crop-specific mask [17,18,19]. For example, Maselli et al. [20] used Normalized Difference Vegetation Index (NDVI) thresholds to isolate the cropland pixels of interest in the Sahelian region. Genovese et al. [18] showed that using the crop class from a medium resolution land cover classification improved coarse resolution yield forecasting accuracies. The authors of [21] used a general cropland mask to extract Moderate Resolution Imaging Spectroradiometer (MODIS) NDVI values for areas with annual crops in the Canadian prairies to forecast barley, canola peas, and spring wheat. Lopez-Lozano et al. [22] used a static coarse resolution percent-cultivated land mask over the EU to extract a crop fAPAR time series for exploring the correlation with official yield statistics. Kowalik et al. [23], using as a reference the arable land fraction from the static Corine Land Cover (CLC, [18]), developed an empirical model to forecast wheat yield over Europe. A general cropland mask can be used to isolate a general cropland signal but does not provide a crop-specific signal, which is preferable for yield forecasting in areas with a diversity of crops.

Due to the widespread practice of crop rotation, the use of a crop-specific mask requires that a season-specific mask be produced to develop the statistical relationship needed for the yield forecast model. For example, Rojas [24] for a study on corn yield in Kenya developed and used crop-specific masks to improve the accuracies of their yield forecasts. Funk and Budde [25] for a study on NDVI-based crop production anomaly estimates in Zimbabwe use cropland mask and NDVI agro-phenological filtering to isolate maize signal and estimate production. Moreover, in order to forecast yields during the growing season, the required crop-specific masks can present a significant logistical challenge [14]. In the USA, Bolton et al. [26] compared the effectiveness of a detailed crop type mask versus a general cropland mask within a study on forecasting corn and soybean yields. Johnson [27] assessed correlations between yields of major US grown crops and a variety of crop-specific time-series of MODIS products including Leaf Area Index, NDVI, and Land Surface Temperature (LST) that were extracted using the USDA NASS Cropland Data Layer (CDL). Guan et al. [28] worked directly at the county scale and used county-level crop yield and area statistics from NASS to select the high-crop counties with >40 total crop area and to derive county crop-yield-based NPP for assessing the US Corn Belt crop yields using a range of satellite data using a partial least squares regression approach. Becker-Reshef et al. [29] developed an empirical but generalized crop yield model based on a single year USDA CDL (2007) crop-specific mask in the US and applied the same approach with a constant crop-specific mask developed for Ukraine. The model was calibrated in Kansas and applied at 0.05° spatial resolution in both Ukraine and the USA. It is based on the relationship between the NDVI value at the peak of the season and the final yield value, corrected for the average purity or the percentage of crop of the 5% purest pixels within the area studied. Later, Franch et al. [30] improved the temporal aspect of the model by including Growing Degree Days (GDD) information. The model was applied satisfactorily (errors lower than 10%) to forecast wheat yield at the national level over the US, Ukraine, and China. It also showed good performance when applied to the AVHRR Long Term Data Record (LTDR) [31].

Even in the case of the USA where the USDA National Agricultural Statistics Service (NASS), a world leader in operational crop type mapping and monitoring, produces annual crop type maps, [32] this information is only publically released after the end of the growing season. Though some promising progress has recently been made in producing early season crop type maps utilizing moderate- and high-resolution data [33,34,35,36], the quality of such maps is usually predicated on the availability of accurate reference (ground truth) data and/or knowledge on crop calendar and phenology, which is often not readily available, especially within a global context. The creation of a single classification model for crop type mapping applicable to multiple years remains a challenge [37] and is an active area of research.

Within this context, this study focused on developing a practical approach to winter wheat yield forecasting applicable at national/sub-national scales in the absence of up-to-date (i.e., within the current growing season) detailed crop-type information. Give the robustness of Becker-Reshef et al. [29] and Franch et al. [30] yield model and our experience with it, we chose to use this approach in this study in order to evaluate the impact of utilizing a single out-of-season crop mask, for yield forecasting. The approach developed was based on spatially aggregating a single, moderate resolution, crop-specific mask from a recent (within 8 years) preceding year, making the assumption that the total area of winter wheat planted varies little (i.e., less than 10%) between successive years. This provides a considerable advantage, as, in general, preceding year masks are more likely to be available than in-season crop-type masks. This is largely due to the fact that detailed crop-type masks are easier to produce once the growing season is complete, as crop phenology (i.e., crop temporal profiles) can be used to differentiate crop types [14,38,39,40,41]. The underlying assumption in this study was that coarsening the scale of the winter wheat mask would result in a constant per-pixel wheat proportion over multiple years, enabling a coherent and comparable extraction of crop-specific vegetation indices for yield estimation over large areas. The approach proposed can have significant implications for operational monitoring at national scales where crop area does not change significantly from one year to the next, as a single, recent, crop-specific mask from a preceding season can be used in combination with coarse spatial resolution NDVI time series for yield forecasting over multiple years.

2. Data and Study Area

2.1. Study Area

The USA was chosen as the study area for this research, as it is one of the main producers and exporters of wheat globally, and because reliable annual winter wheat classifications are available from the USDA NASS Cropland Data Layer (CDL), which can be used to reliably evaluate our approach.

Amongst field crops in the US, wheat ranks third in terms of planted area behind corn and soy. Wheat is produced in almost every state in the United States and winter wheat varieties dominate US production, representing between 70 and 80% of the total. The main class is Hard Red Winter Wheat, which is grown primarily in the Great Plains, with Kansas being the largest producing state. A study by Yan and Roy on field sizes in the USA [42] found that approximately 85% of wheat fields are smaller than 1 km² and that wheat fields are generally larger than the other crop types for sizes equivalent to between 0.5 × 0.5 mile (0.646 km²) and 0.5 × 1 mile (1.295 km²), which are both common US field sizes.

Due to the widespread practice of crop rotations in many parts of the country, the location of planted wheat fields alternates from one year to the next within the majority of counties. Figure 1 shows an example of crop rotation practices in the particular case of Kansas. Though it should be noted that the total planted area remains relatively constant (Figure 2). According to the USDA NASS statistics, the standard deviation of wheat planted area in the US between 2008 and 2017 was 10% [43].

This study relied on two types of available data: (i) Crop type maps from the USDA NASS Cropland Data Layer (CDL) for eight years, 2008 through 2017, and (ii) NASA MODIS 250 m daily NDVI [44,45]. These datasets are summarized below.

2.2. The USDA NASS Cropland Data Layer (CDL)

The CDL is a crop-specific, rasterized, geo-referenced land-cover map produced by NASS [32,46,47]. NASS is the agency responsible for administering the USDA’s U.S. program, collecting and publishing agricultural statistics at the national, state, and county levels. The CDL is produced annually for the whole United States since 2008 (previous years are also available but just cover the main agricultural states) and is released to the public domain after the growing season and release of official county crop estimates. The CDL layers used in this study were produced by NASS using moderate resolution (30–56 m) imagery from the Indian Advanced Wide Field Sensor (AWiFS), Landsat 5, Landsat 7, the NASS June Agricultural Survey (JAS), and Farm Service Agency (FSA) Common Land Unit for ground truth and training data. The major crop types, including winter wheat, within the CDL generally have a classification accuracy “ranging from mid-80% to mid-90%” [48]. In this study, CDL winter wheat binary masks for 2008 through 2017 were extracted and aggregated to a range of increasingly coarser spatial resolutions in order to investigate the relationship between crop rotation and spatial scale, and the feasibility of using spatially aggregated crop masks to forecast yields.

2.3. Vegetation Index Time Series Data

Vegetation indices (VI) such as NDVI have long been recognized for their value for monitoring vegetation changes as inputs to crop yield and production forecasting models, as well as for tracking crop growth and development [8,25,49,50,51,52,53]. In particular, it has been shown that yield-forecasting models are significantly improved when crop masks are used in conjunction with the VI data, as the crop masks isolate the signal of the crop under investigation [14,54]. In this study, the CDL-derived binary winter wheat masks were scaled up to a range of coarser spatial resolutions as percent wheat masks and used to extract NDVI time-series for winter wheat pixels. These time series were used to estimate winter wheat yields using the GLAM Wheat model, refs. [29,30] to evaluate the utility of coarsening resolution for yield forecasting using a single season wheat mask.

Crops grow and develop rapidly during the growing season. Winter wheat is generally planted in autumn and crop emergence and stand establishment typically occurs in late autumn prior to winter. During winter, vegetative growth slows or ceases and the crop’s resistance to cold weather is increased. In the spring, winter wheat resumes rapid vegetative growth followed by reproductive development and reaches maturity in early to mid-summer. As such, frequent observations are critical for effectively monitoring the crop’s rapid development and state, which currently is only achieved by coarse spatial resolution sensors such as MODIS.

Daily surface reflectance at 250 m resolution from MODIS Collection 6 product (M{O,Y}D09GQ) time-series [44] were obtained for the 2008 through 2017 growing seasons. The surface reflectance product is produced by NASA MODIS Advanced Processing System (MODAPS) [55]. Based on the VJB method [56,57], we derived the nadir Bidirectional reflectance distribution function (BRDF)-corrected surface reflectance, which was used to compute the NDVI using the Red (MODIS band 1, 620 nm–670 nm) and near infrared (MODIS band 2, 841 nm–876 nm) bands [58]:

$$\mathrm{NDVI}=\left(\mathrm{NIR}‒\mathrm{RED}\right)/\left(\mathrm{NIR}+\mathrm{RED}\right).$$

(1)

We used daily surface reflectance observations (M{O,Y}D09GQ) with BRDF corrections rather than the standard MODIS 250 m NDVI product (M{O,Y}D13Q1) for two reasons. Primarily, the MODIS NDVI product is only available as a 16-day or monthly composite while the surface reflectance product corrected for BRDF effects provides daily NDVI at a 250 m resolution. The latter therefore improves the likelihood of capturing the maximum NDVI value. Secondly, the MODIS NDVI product is not BRDF corrected; rather, it is based on a series of criteria that selects the cleanest (i.e., least cloud and aerosol contamination) pixel in the compositing window with the best viewing and illumination geometries. The BRDF correction method normalizes the effects of viewing and illumination geometries, providing a more consistent and comparable value between dates.

3. Methods

This study was carried out in four main steps: (1) spatial aggregation of the 2008–2017 winter wheat masks to increasingly coarser spatial resolution percent-wheat masks; (2) analysis of inter-annual percent wheat variability over eight years at a series of spatial resolutions; (3) extraction of winter wheat-specific NDVI temporal profiles; (4) yield estimation under two scenarios for crop mask availability.

3.1. Spatial Aggregation of Winter Wheat Masks

The US CDLs (2008–2017) were used to create eight binary winter wheat masks at the CDL’s native spatial resolution. The native spatial resolution of the CDL for 2008 and 2009 was limited to 56 m due to the use of AWiFS imagery as input. Subsequent versions of the CDL are available at 30m spatial resolution. At the native resolution wheat pixels were assumed to be pure pixels and therefore assigned a value of 100% and non-wheat pixels were set to 0%. Figure 3 shows an example of the 2017 growing season winter wheat mask at 250 m spatial resolution. These masks were then re-projected to match the MODIS NDVI data sinusoidal projection.

Each of the eight annual winter wheat masks was scaled up to a series of seven, increasingly coarser resolution percent wheat masks for a total of 70 masks, where pixel values represent the proportion of wheat within each pixel. This translates into seven masks for each year between 2008 and 2017 at the following spatial resolutions: 250 m, 500 m, 750 m, 1 km, 2 km, 5 km, and 10 km. At each resolution, the per-pixel wheat percent was calculated by averaging the value of the original resolution binary masks (wheat—100, not wheat—0) within the coarser resolution pixel (Equation (2)). This approach is similar to the spatial aggregation approach commonly described in the literature [59,60,61].

For pixel x at coarse resolution r:

$$Pct{W}_{x,r}=\frac{nCDLWheat}{N}\times 100$$

(2)

where $Pct{W}_{x,r}$ is the scaled-up percent-wheat value for a given pixel x, at coarse resolution r; N is the number of 56 m pixels within pixel x; nCDLWheat is the number of wheat pixels derived from the CDL winter wheat mask within pixel x.

3.2. Analysis of Inter-Annual Percent Wheat Variability

At the basis of the developed approach is the assumption that, at a certain spatial resolution, the percent of wheat within a given pixel remains stable from one year to the next (within a limited number of recent years) despite the practice of crop rotations. Implied is that the overall planted area within the region of interest is relatively constant over the specified time frame. In the case of the US, the winter wheat planted area varied by 10% between 2008 and 2017. Maselli et al. [62] made a similar assumption in a study on wheat production estimation in Tuscany, where they assumed wheat covers a constant fraction of winter crop pixels in all provinces and years used in the study. To assess this assumption in the US and determine the spatial resolution at which the per-pixel wheat purity remains constant over a period of 10 growing seasons, the range of the percent wheat values (RPct) for every pixel (x) at all 7 resolutions (r) was computed (Equation (3)).

$$RPct{W}_{x,r}=Max\left(Pct{W}_{x,r}\right)-Min\left(Pct{W}_{x,r}\right)$$

(3)

where $Max(Pct{W}_{x,r})\mathrm{and}Min(Pct{W}_{x,r})$ are respectively the maximum and the minimum percent wheat value for pixel x at resolution r over 10 years.

3.3. Winter Wheat-Specific NDVI Time Series

The next step was to evaluate the effect of using a single season percent winter wheat mask for extracting NDVI temporal profiles over multiple years (2008–2017) for the range of spatial resolutions. Following [29], the average NDVI of the 5% purest pixels (VI_wheat) was adjusted for background noise by subtracting the 5% minimum mean NDVI signal through the time series (VI_background) to create a time series for each county (Equation (4)).

$$\begin{array}{c}ANDVI\left(day,county\right)\\ =NDV{I}_{wheat}\left(day,county\right)\\ -{\displaystyle {\displaystyle \sum }_{year=2008}^{2017}}NDV{I}_{background}\left(year,county\right).\end{array}$$

(4)

In order to select the optimal level of spatial aggregation for yield forecasting in the US (accounting for the dampened signal effect as the resolution was coarsened), the GLAM wheat yield model (fully described in [29,30] was run under 2 scenarios of winter wheat mask availability for the full range of pixel spatial resolutions using the extracted MA_NDVI and percent wheat mask as inputs. The wheat mask availability scenarios are described further in Section 3.5.

3.4. GLAM Wheat Yield Model

The GLAM Wheat model, fully described in [29] and [30] is an empirical, generalized regression model based on the assumption that the yield is positively and linearly correlated to the seasonal maximum ANDVI (ANDVI_peak,year) at the county level and to the purity of the wheat signal. The ANDVI_peak,year is estimated averaging the 1% maximum ANDVI values for each year and county. Following [29], the average of the percent wheat values of the purest 5% winter wheat pixels for each county (Mcpt) is used to derive the percent wheat dependent slopes using Equation (5).

$$S_{Mpct}=9.61+\left(-0.05Mpct\right).$$

(5)

We then applied this slope uniformly to the entire US to obtain winter wheat yields in tons/ha for each county by multiplying the ANDVI_peak,year value by the corresponding slope derived (Equation (6)).

$$ForecastYield=ANDV{I}_{peak,year}{S}_{Mpct.}$$

(6)

In this way, the seasonal maximum ANDVI is weighted by the maximum purity that can be observed in each county. The forecast production was finally obtained by multiplying the yield forecast by the official county area harvested statistics.

$$TotalForecastProduction={\displaystyle \sum }_{AU}ForecastYield\times AreaHarvested.$$

(7)

The model was developed in the state of Kansas (USA) and then tested over Ukraine. Franch et al. [30] enhanced the Becker-Reshef et al. [29] model in order to get an improvement on the timeliness of the production forecast by including the information of the accumulated Growing Degree Days (GDD). The model was applied satisfactorily to MODIS and NCEP/NCAR data over the US, Ukraine, and China, and was proven to be applicable to other sensors i.e., VIIRS [63] and Landsat [35]. Most recently, the yield model was applied to AVHRR historical data from 1982 to 2017, over the U.S. [31]. In all these cases, it was able to estimate winter wheat yields within 10% of official statistics, up to two and a half months prior to harvest. The enhanced winter wheat yield model [29,30] is used in this study to estimate yield in this study.

3.5. Yield Estimation under Two Crop Mask Availability Cases

Ideally, detailed and annually updated winter wheat maps would be available during the growing season at an adequate spatial resolution to extract NDVI time series profiles of the target crop. At this time, in-season crop-specific masks are rarely available for the majority of agricultural regions [26,64]. However, a single or a series of crop type maps are often available or can be more easily generated for the previous season, after the growing season is complete.

This study therefore sought to assess the feasibility of using a previous season’s spatially aggregated mask for yield forecasting and to assess the tradeoff between the “optimal” situation where season-specific winter wheat masks are available during the growing season (Case 1—control case) and the most likely situation where only a single mask from the previous season is available (Case 2).

To accomplish this task, the GLAM wheat yield model was run for the two cases at all spatial resolutions. Case 1 was treated as the “control” case representing “optimal” data availability. As such, Case 2 results were assessed relative to the control case. The estimated yields were also assessed against the official NASS reported yields and the root mean square error (RMSE Equation (8)) and corresponding percent error were computed for each case and at each spatial resolution. At the time of this study the NASS winter wheat masks for the whole US were available for 2008 through 2017 so Case 1 was evaluated on these years.

$$RMSE=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(Yiel{d}_f-Yiel{d}_o\right)}^2}{n}}$$

(8)

where n is number of years, Yield_f is the forecast yield, and Yield_O is the observed yield.

In Case 2 (which represents instances when only the previous’ season winter wheat mask is available), the previous season’s mask was used to extract the NDVI signal. In other words, the 2008 mask was used to extract the 2009 growing season NDVI signal and so on. As for Case 1, the RMSE was computed for all Case 2 model runs.

4. Results

4.1. Analysis of Inter-Annual Percent Wheat Variability

A total of 70 percent-wheat masks were created, 7 masks for each year between 2008 and 2017 at the following spatial resolutions: 250 m, 750 m, 500 m, 1 km, 2 km, 5 km, and 10 km. These masks were then used to assess the variability of the per pixel percent wheat at each spatial resolution.

Figure 4 presents a series of images focusing on Western Kansas depicting the RPctW values over 10 years in an area with high crop rotation rates. At finer resolutions (i.e., 250 m), RPctW values were close to 100, indicating that in some years no winter wheat was planted in that location and in other years the pixel was purely winter wheat (Figure 4 shown in red-orange). As pixel size was increased, the variability of percent wheat within a given pixel decreased and stabilized.

To quantify and characterize this percent-wheat “stabilization” effect over the entire country, the median RPctW (MRPct) was computed for all considered counties (Equation (9)) at each of the spatial resolutions and is presented in Figure 5. These results indicate that the percent wheat inter-annual variability drops below 10% after data are aggregated to 5 km for the majority of counties in the US.

$$MRPc{t}_{c,r}=Median\left(RPct{W}_{x,r}\right)|x\in c$$

(9)

where $MRPc{t}_{c,r}$ is the Median Range Percent Wheat for county c at resolution r.

4.2. Winter Wheat-Specific NDVI Time Series

Figure 6a,b shows the results of extracting winter wheat-specific NDVI time series for two contrasting types of counties. Figure 6a for Decatur County, Kansas, where crop rotation is common, and Figure 6b for Harper County, Kansas, where winter wheat monoculture is widespread. Each time-series graph depicts the MA_NDVI signal from 2008 to 2017 extracted using the 10 annual crop masks produced in this step, at all spatial resolutions. In the case of Decatur County, it is evident that a single, season-specific crop mask cannot be used to extract the MA_NDVI at the higher resolutions due to the effect of crop rotations. However, as the resolution is coarsened the extracted signal from each of the crop masks converges, indicating that at coarser resolutions a single mask can be used over multiple years to extract a consistent signal. On the other hand, in Harper County almost all the masks, even at the finer resolutions, can be used to extract a consistent signal since winter wheat fields are not heavily rotated between years.

4.3. Yield Estimation under Two Crop Mask Availability Cases

The regression-based winter wheat yield model was run for the cases described. The results for Cases 1 and 2 for the US are presented in Figure 7. Case 1 (red in Figure 7), as expected, had the overall lowest errors (5–6%) with a stable model error across all resolutions given the statistical significance of the results represented by the error bars shown in Figure 7. Case 2 performed poorly relative to Case 1 at the finer resolutions due to the effect of crop rotations, though performed remarkably well as the spatial resolution was coarsened. The error was highest (25–30%) at the finest resolutions, and rapidly declined with coarser spatial resolution with errors equivalent to Case 1 at the 5 and 10 km resolutions.

Figure 8 shows an example of the percent error evolution analogously to Figure 7 but for two contrasting states Kansas and Washington. Kansas shows a lower difference between Case 1 and Case 2 while Washington state shows a larger difference between both cases that can increase up to 50% for the finest spatial resolutions. The reason for such difference in Washington (compared to Kansas) can be explained by the higher overall rate of crop rotation in this state. In both cases, the results provide the same conclusion that the national level plot, with errors from Case 1 and Case 2 equivalent at coarser spatial resolutions.

These results are quite promising as they indicate that, in the absence of a timely seasonal crop-type mask, a single percent wheat mask, derived from a recent season’s crop-type mask, can be used in combination with coarse resolution NDVI to forecast yields prior to harvest. The tradeoff of using a single mask in terms of percent error relative to the Case 1 minimum error (5.4–5.9%) is 0.2–0.7% in Case 2.

5. Discussion and Conclusions

Crop yield forecasting is hampered due to limited availability of detailed crop type masks prior to the end of the growing season. The objective of this study was to develop and evaluate a practical approach to winter wheat yield forecasting that utilizes readily accessible data to forecast yields at sub-national and national scales with minimized tradeoff in accuracy in cases where availability of crop type masks are limited.

To this end, an approach utilizing a single static crop-type mask derived for a preceding season/s, with freely available coarse resolution NDVI time series was evaluated in the US. The approach is based on spatially aggregating moderate resolution (30–56 m) winter wheat masks, to a spatial resolution at which the per-pixel percent wheat remains constant between years enabling extracting a coherent crop-specific NDVI signal from pixels covered with varying proportions of winter wheat. These data were used to run a yield model that utilizes winter wheat percent and maximum seasonal NDVI to forecast yields.

This study found that the best yield estimates were attained when moderate resolution winter wheat masks were available for every season analyzed. However, the tradeoff between using a coarse-resolution, static mask derived from a single preceding season was surprisingly low (0.2–0.7% for 10 and 5 km resolutions, respectively) when compared to the (more optimal) mask representing the current growing season. These results from the US suggest that, in the absence of within-season, detailed crop-specific masks, yield can be forecast at the national level with a small tradeoff in accuracy, using readily available coarse resolution data and a preceding season’s percent wheat mask.

It should be noted that the model provides errors lower than 6% for most spatial resolutions and using in-season crop masks. However, the accuracy of the model decreases for high spatial resolutions, reaching maximum error at 500 m. This is likely due to several factors including the accuracy of the crop masks, artifacts of 250 m gridding such as pixel growth, shift, and distortion (i.e., it has been shown that at the 250 m resolution there can be up to a 0.5 pixel shift [64,65,66,67]). In addition, the simple aggregation approach implemented in this study does not account for the MODIS spatial response, whereby off-nadir sensor acquisitions, i.e., those at increasingly higher view angles, will have a ground sample size greater than acquisitions at the nadir viewing angle (the so-called “bow tie effect”). As pointed out by Duveiller et al. [68], this can impact the percent winter wheat characterization in the aggregated crop type masks. Thus, this topic should be further explored.

While one would assume that, for agricultural remote sensing, it is best to work at the finest spatial resolution so that individual fields can be resolved, this study provides a counter-intuitive result. The results suggest that coarse spatial resolution data with high temporal frequency, combined with lower temporal resolution moderate spatial resolution data, can offer a viable alternative in the absence of finer spatial data and up-to-date crop-specific masks, with a small tradeoff in accuracy in this case. The noise in the crop signal introduced by year-to-year variability due to crop rotations is reduced by monitoring crops at a coarser resolution. This is in keeping with results presented by Markham and Townshend [69], Malingreau and Belward [70], Nelson et al. [61], and more recently Leroux et al. [71], albeit at finer spatial resolutions and for different land cover types.

This study presents a viable and efficient method for isolating a coherent, crop-specific signal over multiple growing seasons using a static, fine resolution crop mask aggregated to a coarse resolution. A fundamental condition for this approach is that there are no major shifts in planted area between years due to such pressures as biofuel demand or market incentives. Nevertheless, according to USDA national level statistics, this condition seems to hold true, at least at the national level, for many of the main wheat exporting countries where the variation in planted area over the past 10 years is below 10% [72]. In the case of the US, where 10 annual national winter wheat masks were available from the USDA NASS, the analysis suggests that over a period of at least 10 years this assumption holds true and that per-pixel percent wheat variability in the US stabilized (below 10%) at approximately 5 km pixel resolution.

The main advantage of the results presented for crop monitoring and specifically for yield forecasting models that rely on remotely sensed inputs is that a single, coarse resolution percent wheat mask is sufficient to extract a consistent winter wheat temporal profile over multiple years. Therefore, a single mask from a previous year would likely suffice. Second, this method helps to mitigate the challenge of generating a timely and accurate crop type mask prior to the end of the growing season in order to forecast yields.

We highlight that the focus of this approach is for forecasting yield over large regions (i.e., at the state or national levels) rather than for precision agriculture purposes or for seasonal crop area estimates. The approach of scaling up a representative seasonal crop-specific mask or a temporally aggregated mask as a percent crop mask is primarily advantageous in regions where crop rotations are common practice. In regions where crop monocultures are prevalent, a single year, moderate resolution crop mask would be as effective. As a coarse resolution crop type mask results in mixed pixels with various crop percentages, a yield model using this approach must explicitly account for the crop percentage as demonstrated in [8] or in Rasmussen et al. [73] as a way of isolating the crop-specific signal. Furthermore, this study only explored such an approach for winter wheat, which has the advantage of greening up prior to summer and spring crops, so further research is needed in order to assess its applicability to other major crop types.

The recent advances in the readily available moderate resolution (20–30 m) data every 3–5 days globally (owing to Sentinels and Landsat), alongside the proliferation in high performance computing and cloud systems, and advances in big data analytics will undoubtedly move the agricultural monitoring community forward, including towards the goal of openly available in-season, moderate-resolution crop type masks. The development of operational capabilities for producing such timely crop type masks, is a high priority for the agricultural monitoring community, as highlighted by the G20 GEOGLAM initiative priorities. At the moment, although significant progress has been made in this direction particularly for the major production regions of the world, (i.e., the Sen2Agri project [74]), these capabilities are still mostly in the applied R&D stages. As such, until such masks are reliably available, the proposed solution in this research seems promising in many of the major production regions (as well as food-insecure ones).

Due to the limited number of seasonal wheat masks that were available to carry out the analysis, it is recognized that the results may be inconclusive and that further assessment should therefore be carried out prior to drawing decisive conclusions. Nevertheless, these results are encouraging and could have promising implications for developing practical, operationally viable methods for crop forecasting at national levels using EO in cases where year-to-year variability of crop area is limited.