Next Article in Journal
Super-Resolved Multiple Scatterers Detection in SAR Tomography Based on Compressive Sensing Generalized Likelihood Ratio Test (CS-GLRT)
Previous Article in Journal
Land-Cover Classification of Coastal Wetlands Using the RF Algorithm for Worldview-2 and Landsat 8 Images
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Evaluation of the Uncertainty in Satellite-Based Crop State Variable Retrievals Due to Site and Growth Stage Specific Factors and Their Potential in Coupling with Crop Growth Models

1
Department of Electrical Engineering, City College of New York, 160 Convent Ave., New York, NY 10031, USA
2
Department of Geography, University of Wisconsin-Madison, 550 N. Park St., Madison, WI 53706, USA
3
Nelson Institute Center for Sustainability and the Global Environment, University of Wisconsin-Madison, 1710 University Avenue, Madison, WI 53726, USA
4
Department of Forest and Wildlife Ecology, University of Wisconsin-Madison, 1630 Linden Drive, Madison, WI 53706, USA
5
CNR-Institute of Mediterranean Forest and Agricultural Systems, 85 Via Patacca, 80040-Ercolano (Napoli), Italy
6
Department of Electrical and Computer Engineering Technology, Farmingdale State College, 2350 Broadhollow Road, Farmingdale, NY 11735-1021, USA
*
Author to whom correspondence should be addressed.
Remote Sens. 2019, 11(16), 1928; https://doi.org/10.3390/rs11161928
Submission received: 22 July 2019 / Revised: 12 August 2019 / Accepted: 15 August 2019 / Published: 17 August 2019
(This article belongs to the Section Remote Sensing in Agriculture and Vegetation)

Abstract

:
Coupling crop growth models and remote sensing provides the potential to improve our understanding of the genotype x environment x management (G × E × M) variability of crop growth on a global scale. Unfortunately, the uncertainty in the relationship between the satellite measurements and the crop state variables across different sites and growth stages makes it difficult to perform the coupling. In this study, we evaluate the effects of this uncertainty with MODIS data at the Mead, Nebraska Ameriflux sites (US-Ne1, US-Ne2, and US-Ne3) and accurate, collocated Hybrid-Maize (HM) simulations of leaf area index (LAI) and canopy light use efficiency (LUECanopy). The simulations are used to both explore the sensitivity of the satellite-estimated genotype × management (G × M) parameters to the satellite retrieval regression coefficients and to quantify the amount of uncertainty attributable to site and growth stage specific factors. Additional ground-truth datasets of LAI and LUECanopy are used to validate the analysis. The results show that uncertainty in the LAI/satellite measurement regression coefficients lead to large uncertainty in the G × M parameters retrievable from satellites. In addition to traditional leave-one-site-out regression analysis, the regression coefficient uncertainty is assessed by evaluating the retrieval performance of the temporal change in LAI and LUECanopy. The weekly change in LAI is shown to be retrievable with a correlation coefficient absolute value (|r|) of 0.70 and root-mean square error (RMSE) value of 0.4, which is significantly better than the performance expected if the uncertainty was caused by random error rather than secondary effects caused by site and growth stage specific factors (an expected |r| value of 0.36 and RMSE value of 1.46 assuming random error). As a result, this study highlights the importance of accounting for site and growth stage specific factors in remote sensing retrievals for future work developing methods coupling remote sensing with crop growth models.

Graphical Abstract

1. Introduction

1.1. Background

Mechanistic crop growth models temporally predict the growth of crops as a function of genotype x environment x management (G × E × M) factors [1]. By mechanistically modeling the effects of G × E × M factors and their interactions, crop growth models are able to integrate information about the properties of the seed (genotype), the decisions farmers make both at planting and within the season (management), and the variability in the weather and soil (environment). Examples of these factors in each category of G × E × M are shown in Table 1 [2,3]. In addition to these G × E × M factors, biotic stresses—such as weeds, pests, and diseases—can further limit the growth of crops and these factors are difficult to model, although some recent advances have been made [4]. Nevertheless, in highly developed cropping systems, such as the US corn belt, fields tend to be well-managed and the reduction in yield caused by unmodeled factors, such as biotic stresses, is generally 20% or less [5,6]. As a result, mechanistic crop growth model simulations are able to provide valuable information with relatively strong predictive performance in highly developed cropping systems [6,7].
Assimilation of remote sensing data into crop growth models can be used to reduce the uncertainty in the G × E × M factors (which control crop growth) via calibration [8,9,10,11]. In the calibration approach to remote sensing data assimilation, the model parameters and G × E × M factors affecting crop growth are adjusted by reinitialization until the crop growth model output agrees with the remote sensing observation (as opposed to the updating or forcing approaches where the crop model state variables are themselves directly altered) [9]. However, uncertainty in the remote sensing retrievals of crop state variables, such as leaf area index (LAI), leads to significant challenges [9] in the calibration and determination of the G × E × M factors. This is because the interactions of G × E × M factors in crop growth models are highly non-linear and careful application of inversion techniques is required to determine input parameters from observations [12,13]. As a result, even small uncertainties in the remote sensing retrievals can propagate into significant errors in the G × E × M factors determined by calibration [14]. Therefore, calibration of crop models with remote sensing data is primarily used to analyze output variables, such as yields and biomass, discarding the G × E × M factors determined by calibration as an intermediate step [8,15,16,17,18].
Nevertheless, improved understanding of the G × E × M factor variability can greatly improve our ability to use crop growth models at the regional scale [6,19,20] to predict into the future and answer questions about climate change [21], agricultural policies [22,23], and yield gaps [24]. At the regional scale, G × E × M parameter uncertainty is even more significant due to a lack of calibration data as compared to the field-scale [1,25]. Thus, constraints from measurements other than yield are vital for further reduction in the uncertainty [25] at this scale. Illustrating this point, ref. [25] found that the majority of the uncertainty in LAI simulations for regional simulations of Indian groundnut was parametric uncertainty, indicating the potential of reductions in the uncertainties of satellite retrievals (such as those of LAI) to significantly improve our understanding of G × E × M variability in calibration of regional crop models [26].
The crop state variable retrieval uncertainty is in a large part caused by the variability in secondary factors [27,28,29,30,31,32] that influence the remote sensing measurements, such as cultivar type, soil background, canopy structure, and inherent leaf properties; most of these secondary factors are strongly dependent on site and growth stage [33,34,35,36]. Physical canopy radiative transfer models, such as PROSAIL [37], provide a theoretical model to understand the effect of the secondary factors by forward modeling the top-of-canopy reflectance spectrum from variables describing the soil background, canopy structure, and leaf properties [9]. However, inversion of canopy radiative transfer models is ill-posed [38] and requires the use of a priori constraints to perform the retrievals [39,40]. While temporal [40,41,42] and spatial [40,43] constraints can be used to address the ill-posedness of the retrieval, they are not sufficiently powerful to remove the uncertainty. As a result, assumptions must be made about the canopy structure and leaf properties [40]. Unfortunately, although both canopy structure and leaf properties have a significant effect on the uncertainty of the retrieval [32], it is difficult to constrain them beyond finding appropriate ranges for the values based on land cover [44] and selecting vegetation indices with greater sensitivity to the variable of interest [32,45,46]. However, even though the full spectral modeling can optimize the best choice of vegetation indices for given applications, using vegetation indices in the retrievals directly still results in valuable spectral information being lost, undercutting the benefits of the possibility of using the full spectral information available with canopy radiative transfer models in the retrieval itself [47] as full-spectrum methods have shown good results in the literature [48,49].
However, because of the lack of information available to remove the uncertainty about secondary factors, physical radiative transfer approaches have not dominated over empirical approaches, although these often do not use the full spectral information available from the sensor and lack a theoretical basis to control secondary factors [27,28,29]. The empirical algorithms overcome these issues by directly using training data to learn to use the “subtle spectral features to reduce undesired effects” [47] that make vegetation retrievals difficult. In addition, in some cases, empirical methods are also able to improve the retrievals with auxiliary information [29,50,51].
In empirical approaches, the uncertainty caused by the variability in secondary factors manifests as the “one place, one time, one equation” issue [27] where regressions between the satellite measurements and the crop state variables trained on one set of sites and times do not generalize well to another set of sites and times [27,28]. The issue occurs because most empirical studies develop a global regression relating the satellite measurements to the crop state variables which does not account for the spatiotemporal variability in the secondary factors, although some studies have attempted to use the secondary factors to improve the retrieval [29,50,51]. Specifically, refs. [50,51] find that developing separate regression models for different growth stages provides the best results, while [29] finds that including cultivar, planting pattern, and growth stage in the model could improve the performance of the retrievals. While the secondary factors in [29,50,51] do not correspond to the secondary factors in physical radiative transfer models such as PROSAIL, their indirect connection to the leaf and canopy parameters used by PROSAIL [33,34,35,36] allows them to reduce the uncertainty caused by the secondary effects. Nevertheless, the work on including secondary effects is quite limited and hampered by lack of available data [28] to span the large spatiotemporal variability in these secondary factors, calling for new approaches to address this issue.
In order to address the uncertainty caused by secondary factors, it is necessary to obtain data that covers the extent of their spatiotemporal variability. Crop growth models provide one possible avenue to obtain information on the secondary factor leaf and soil properties. The use of crop growth models to obtain information about the secondary factors has been best explored in coupling studies [52,53,54,55], where remote sensing data is assimilated into a combined model consisting of a crop growth model, a canopy radiative transfer model, and formalisms linking the outputs of the growth model with the inputs of the radiative transfer model. These studies [52,53,54,55] have been successful in coupling several variables from the crop growth models, such as LAI, leaf structure parameter, water content, dry matter content, total chlorophyll content, and relative soil dryness. The variables coupled in addition to LAI are secondary factors that affect LAI retrieval [32] and the coupling can be understood to provide constraints on these secondary factors from the biological mechanics of growth and its interaction with the weather/soil environment. In addition, if available, any genetic (cultivar choice) or management information inputted into the crop model can provide additional constraints on the secondary factors [56]. Unfortunately, it is difficult to use crop growth models to gain information about these secondary parameters at a regional scale as information about G × M parameters is limited at this scale [57]. As a result, regional crop growth model simulations are generally validated only against crop yields and phenological dates [6,20,58,59,60] and consequently may have significant uncertainty in their prediction of in-season state variables (many of which are secondary factors in LAI retrieval) [61]. In contrast, field-scale crop growth model simulations have been validated in much more detail with respect to in-season state variables. For example, several studies [2,62,63,64,65] evaluate their performance in predicting LAI, canopy cover, biomass, soil moisture, soil nitrogen, plant nitrogen, evapotranspiration, and phenology as well as yield. The crop model’s stronger performance at field-scale in predicting both the yield and individual within-season process can be attributed to the availability of significantly more accurate agromanagement information, and to a lesser extent to more accurate soil and weather data, at this scale [66]. Thus, incorporating field-scale crop growth modeling of secondary parameters in training and testing agricultural satellite retrieval algorithms [67] can potentially provide for significant advances in addressing the uncertainty caused by site and growth stage specific secondary factors.

1.2. Overview

In this study, we seek to show that the difficulties in using remote sensing to determine the G × E × M factors affecting crop growth are strongly connected to variability in the relationship of satellite measurments and crop state variables and that the variability in the relationship is in a large part caused by site and growth stage specific factors. In order to achieve these objectives, this study uses field-scale crop growth model simulations powered by accurate agromanagement information and collocated with satellite data at the Mead, Nebraska Ameriflux sites, supplemented by ground-truth data from additional sites for validation. Crop growth model simulations are used from only the Mead, Nebraska Ameriflux sites because geolocated agromanagement information, vital [66] to strong simulation performance, is difficult to collect, partially due to farmer concerns about data privacy [68], limiting available information about commercial-sized plots. The availability of collocated crop growth model simulations allows us to (a) analyze the sensitivity of the genotype x management (G × M) factors retrieval by the satellite to variability in the relationship of satellite measurments and crop state variables and (b) use time-series analysis to analyze the uncertainty caused by this variability. Furthermore, the collocated crop growth model simulations are used to demonstrate the possibility of training and testing agricultural remote sensing algorithms with farmer-collected agromanagement data across a wide range of spatiotemporal variability, following the concept we introduced in [67] at the regional scale. Specifically, as in [67], the crop growth model simulations based on the provided data can be used to train and test remote sensing retrieval algorithms and, with sufficient farmer participation, a large swath of the spatiotemporal variability of the secondary factors affecting the retrievals can be covered. This dataset would allow further research to find methods to optimally use available weather, soil, and remote sensing data to create algorithms to map the regional-scale variability in G × E × M. As a result, by using crop growth model simulations at a fixed number of sites where the G × M parameters are known, a remote sensing retrieval algorithm could be trained to map G × M parameters where they are unknown and where no high quality collocated crop growth model simulations are available.

2. Materials and Methods

2.1. Data

In this study, we rely on two ground-truth maize datasets, which we term FLUX and LAIGROUND. The data sources are summarized in Table 2.
The FLUX dataset consists of CO2 flux tower measurements of gross primary productivity (GPP) and incoming solar radiation (SRAD) time series in maize fields. The eddy-covariance technique determines the CO2 flux, which is termed the net ecosystem exchange (NEE), from the covariance of the vertical wind velocity and CO2 flux, sampled by the tower at 10–20 Hz and averaged to 30–60 minute periods [69]. The height of the flux tower is selected to have an appropriate footprint covering the field being studied by the tower. The ecosystem respiration is removed from the NEE to obtain the amount of carbon captured by the producers in the field (GPP) by a partitioning algorithm. In this study, the GPP is either obtained from the nighttime-partitioned product provided by FLUXNET2015 [70] or the site principal investigators (PIs), or calculated from NEE using the nighttime-based partitioning algorithm of [71] implemented in [72]. In addition, ground-truth LAI that was measured at sites on some days of the season and the planting and harvest dates were obtained.
The LAIGROUND dataset consist of ground-truth LAI measurements of maize obtained during various campaigns with different measurement technique (Destructive, LAI2000, AccuPAR, Hemispheric Photography) compiled by [27]. Destructive measurements of LAI rely on physically sampling leaves in predefined areas in the field and measuring them in a laboratory to estimate the LAI in the field. In contrast, the LAI2000, AccuPAR, and Hemispheric Photography techniques use ground-based optical measurements made by researchers in the field on sampling campaign days, along with physics and image-processing based techniques, to estimate the LAI. Further details on all the different measurement techniques can be found in [73]. Each site in this dataset represents a different measurement campaign and some consist of LAI measurements on a single day in neighboring plots, some consist of LAI measurements in different fields (sometimes many kilometers apart), and some consist of multitemporal measurements in the same field/plot. Two of the sites are taken at CO2 eddy-covariance tower sites in the FLUX dataset (Italy and Mead) and the analysis conducted in this study takes care to ensure these are treated as the same sites across datasets when any site-based cross-validation-type analysis is conducted. Following [27], LAI measurments greater than 6 and less than 0.1 are excluded from the LAIGROUND dataset as they are beyond the prediction power of vegitation indicies.
In addition to the ground data in Table 2, we also use solar-reflective satellite data collocated with the ground data. Data from the Thematic Mapper (TM) sensor was used from LANDSAT 5, while data from the Enhanced Thematic Mapper Plus (ETM+) sensor was used from LANDSAT 7. The LANDSAT satellites used for each site depend upon which LANDSAT satellites were active when the site’s data was collected; LANDSAT 5 was active from March 1984–January 2013, while LANDSAT 7 was active from April 1999 to present (ca. August 2019). Data from both satellites was used at sites where data was collected when both satellites were active. For the LAIGROUND dataset, the plots tend to be small and we consequently use 30-m atmospherically-corrected LEDAPS surface reflectance data from LANDSAT 5 and 7 obtained from Google Earth Engine via the GEEXTRACT python tool within 5 m of the plot coordinates. For the FLUX dataset, the plots tend to be production-sized fields and we obtain the average LANDSAT LEDAPS [74] surface reflectance within a 100-m radius of the plot coordinates. In addition, because the LANDSAT temporal resolution is quite low, we obtain MODIS MCD43A4 BRDF-corrected nadir surface reflectance [75] at daily time steps (based on a weighted window of 16 days of measurements) at 500 meters for the FLUX sites, allowing for temporal analysis of the retrieval performance. MODIS data was available for the entire study period for the FLUX sites.

2.2. Hybrid-Maize (HM) Simulations

Simulations from the Mead, Nebraska Ameriflux sites performed by [90] with the Hybrid-Maize (HM) crop growth model are used in this study. The simulations in [90] are based on accurate weather, soil, and agromanagement inputs at the sites and were publicly released [91]. The agromanagement inputs that were recorded at the sites and included in the simulations are planting date, cultivar maturity, plant density, and irrigation. The simulations were validated by [90] with respect to yield, crop respiration, soil respiration, and ecosystem respiration; they are further validated by us in Section 3.1 with respect to LAI and canopy light use efficiency (LUECanopy).

2.3. Methods

In this subsection, we discuss the methods we use to evaluate the influence of site and growth stage specific secondary factors on the relationship between crop state variables and satellite measurments and the retrievability of G × M factors from satellite data. We focus on LAI and GPP in this study because these variables are some of the most commonly retrieved from remote sensing [92]. GPP also serves as a good complement to LAI because, unlike LAI, it is measured on a daily time scale at CO2 eddy-covariance tower stations. Thus, it can be used to provide validation of the temporal analysis performed on crop growth model simulations of LAI. In addition, it should be noted that, as in [67], the methods in this paper can be applied to crop growth model simulated variables whose time series are more difficult to measure than LAI and GPP, providing a basis to analyze performance over a wide range of crop state variables.
As daily GPP strongly depends on the daily SRAD, studies analyzing satellite-derived GPP must account for the strong temporal variability of SRAD when performing retrievals; this is because the variability in SRAD can mask the much smaller variability component in GPP caused by changes in the leaves, plants, and canopy structure [93]. A common technique to do so is correlating the product of the remote sensing measurement and SRAD with daily GPP, as opposed to the remote sensing measurement itself [93]. To achieve a result identical to [93], we analyze the canopy light use efficiency (LUECanopy) in place of the GPP, which we define as
L U E C a n o p y = G P P S R A D ,
As the definitions of various light use efficiencies are not standardized in the literature, we need to clarify that LUECanopy is essentially equivalent to LUEInc in [94], except that incident photosynthetically active radiation (PARinc) is used in place of SRAD. In addition, we wish to note that for the purposes of this study, the criticism of LUEInc in [94] does not apply because our goal in calculating LUECanopy is simply to remove the influence of SRAD and not any plant-based process.

2.3.1. Evaluation of HM Simulations

First, in order to use the HM simulations to evaluate the retrievals, we expand upon the validation performed by [90] to include LAI and LUECanopy. To do so, the modeled and measured values are scatter plotted against each other and the coefficient of determination (R2) to the best-fit line and the root mean square error (RMSE) between the modeled and measured data are calculated. In order to facilitate comparison between the modeling performance of LAI versus LUECanopy, only dates on which both LAI and LUECanopy measurements were available were included in the analysis to ensure that the distribution of crop growth stage did not vary between scatterplots or performance metrics (R2 and RMSE).
In addition, because daily LUECanopy measurements were available, a separate analysis of the performance of the LUECanopy values and the change in LUECanopy is made. The change in LUECanopy is defined as
Δ L U E C a n o p y [ t ] = L U E C a n o p y [ t + Δ 1 ] L U E C a n o p y [ t Δ + 1 ] ,
where Δ is in days and termed the Δ window. ΔLUECanopy is more sensitive to environmental-induced changes than the LUECanopy value itself and the performance in modeling it thus provides additional information on the strengths and limitations of the model.
Furthermore, because of high frequency variability in LUECanopy, the time series modeling performance is analyzed at various levels of smoothing. The smoothing is performed by a moving average filter which is defined as
L U E C a n o p y ¯ [ t ] = 1 2 N 1 i = N + 1 N 1 L U E C a n o p y [ t + i ] ,
where N is in days and termed the smoothing window.

2.3.2. Regression-Based LAI and LUECanopy Retrieval

Second, we train a regression of LANDSAT measurements to LAI and LUECanopy with the LAIGROUND and FLUX datasets. Specifically, we determine the regression coefficients in
L A I = a E V I 2 + b ,
L U E C a n o p y = c E V I 2 + d ,
where EVI2 is the Enhanced Vegetation Index 2 [27] and is defined as
E V I 2 = 2.5 N I R R e d 1 + N I R + 2.4 R e d ,
and NIR is the surface reflectance in the near-infrared band, while Red is the surface reflectance in the red band. The NIR is designated as Band 4 (0.77–0.90 µm) on Landsat 5 and 7, while the Red is designated as Band 3 (0.63–0.69 µm). The coefficients are determined with leave-one-site-out cross-validation by calculating the coefficients on all sites except the one being evaluated. The RMSE performance is then assessed using the coefficients determined from all the other sites and the procedure is repeated for each site. In addition, confidence intervals for the coefficients are determined by bootstrapping. Specifically, for each left-out site, regression coefficients are determined for 1000 random subsets of the remaining sites with the probability of inclusion of a point in any individual random subset equaling 50%. The 5th and 95th percentiles for the regression coefficients of these subset realizations are used as the estimated lower and upper bound of the leave-one-out regression coefficients for the site.
The LAIGROUND and FLUX datasets are analyzed separately for this procedure. The nearest cloud-free LANDSAT measurement within 15 days of the ground measurement is used to analyze the LAIGROUND dataset for consistency with [27], while the average cloud-free LANDSAT measurement within 10 days of the ground measurement is used for the analysis of the FLUX dataset.

2.3.3. Satellite Retrieval and Crop Growth Model Sensitivity Analysis

Third, we analyze the sensitivity of the crop growth model to its G × M inputs and analyze how uncertainty in the satellite retrieval of LAI propagates to the uncertainty in estimation of its G × M inputs. Specifically, we perform new Hybrid-Maize simulations based on the inputs used in [90], varying the planting density, the planting date, and the seed’s growing degree days to maturity from their actual values, and observe the error in the modeled LAI with respect to the measured LAI for the modified simulations. As the emergence date is directly input into the simulations in [90], a preliminary set of Hybrid-Maize simulations is used to determine the appropriate planting date in Hybrid-Maize for the observed emergence date and then this planting date is varied in the sensitivity analysis. This method of determining the planting date to be varied is used in place of the actual planting date to remove the uncertainty caused by modeling the planting to emergence time (as in [90]).
Comparison of the modeled LAI is performed with both the actual measured ground-truth LAI and the measured LAI retrieved from the MODIS measurements. To visualize the effect of the uncertainty in the regression coefficients, the error is shown for a range of regression coefficients determined from the confidence intervals obtained by bootstrapping in the previous subsection. Specifically, the slope of the regression is linearly varied from its minimum lower bound to its maximum upper bound while the intercept of the regression is simultaneously varied from its maximum upper bound to its minimum lower bound. As a large value for the intercept compensates for a lower value in the slope and vice versa, this method generates a realistic space within which to analyze the variation of the regression coefficients.

2.3.4. Evaluation of Uncertainty of LAI and LUECanopy Retrievals Due to Site and Growth Stage Specific Factors with Temporal Analysis

Fourth, we assess the uncertainty of LAI and LUECanopy retrievals with temporal analysis due to site and growth stage specific factors. Due to the “one place, one time, one equation” concept [27], different regression equations should be used to retrieve the LAI and LUECanopy at different sites and growth stages (different times). Furthermore, data from different years may also appear to require different regression equations because the interannual difference in weather and agromanagement is very significant [13] and can cause large differences in secondary factors. Therefore, different years can also be considered different sites for the purposes of this analysis. In order to separate uncertainty caused by site and growth stage specific factors from other types of uncertainty, we use temporal analysis and focus on the retrieval of the temporal change in LAI and LUECanopy. Errors caused by site and growth stage specific factors should be strongly positively correlated at the same place and nearby times; as a result, errors should partially cancel out when retrieving the temporal change as opposed to the actual values themselves. Thus, in order to assess the extent of the uncertainty caused by site and growth stage specific factors, the retrieval error of the change in LAI and LUECanopy is compared to the theoretical error of the change in LAI and LUECanopy assuming temporal independence of error.
To perform the temporal uncertainty analysis for LAI, we use the LAIGROUND dataset as the baseline retrieval and apply the LANDSAT-trained leave-one-site-out regression coefficients from Equation (4) to the MODIS MCD43A4 BDRF-adjusted daily surface reflectance time series to obtain retrievals of LAI with daily resolution. The NIR band is designated as Band 2 on MODIS (0.84–0.88 µm), while the Red band is designated as Band 1 on MODIS (0.62–0.67 µm). The training of the LAI retrieval algorithm is performed on the LAIGROUND dataset with LANDSAT measurements for two reasons:
  • Using the LAIGROUND dataset with LANDSAT imagery better allows for the use of exact point measurements in fields and is thus less likely to be subject to uncertainty in training due to the inhomogeneity of LAI in the field, which can be significant [95].
  • Training on high-resolution LANDSAT imagery as opposed to moderate-resolution MODIS imagery is preferable due to the significance of the mixed-pixel effect and neighboring pixels of other land types (including other crops) [95,96].
In addition, a scaling effect correction algorithm is not used to correct for the uncertainty in applying a regression trained on LANDSAT data to MODIS data as these algorithms generally require a priori information on the subpixel contents of the moderate resolution MODIS pixels [95,96] which is not readily available. For this reason, training on MODIS pixels would likely not provide a benefit with respect to the uncertainty as it is likely that the bias caused by LAI inhomogeneity and the mixed pixel effect varies strongly from site to site [95,96].
With these daily LAI retrievals from MODIS measurements, we calculated the change in LAI as
Δ L A I [ t ] = L A I [ t + Δ 1 ] L A I [ t Δ + 1 ] ,
where Δ is in days and termed the Δ window.
The MODIS-retrieved ΔLAI is compared to the crop growth model predicted ΔLAI using the correlation coefficient absolute value (|r|) and RMSE. These metrics are compared to the theoretical |r| and RMSE if the error of retrieved LAI [t + Δ − 1] and LAI [t − Δ + 1] were independent with a RMSE equivalent to the leave-one-site-out RMSE calculated in Section 2.3.2. In this case, the theoretical RMSE and |r| can be calculated as
R M S E ( Δ L A I [ t ] ) T h e o r = R M S E ( L A I [ t + Δ 1 ] L A I [ t Δ + 1 ] ) = 2 R M S E ( L A I [ t ] ) ,
| r   ( Δ L A I [ t ] ) T h e o r | = | c o v ( Δ L A I a c t u a l + e Δ L A I , Δ L A I a c t u a l ) v a r ( Δ L A I a c t u a l + e Δ L A I ) v a r ( Δ L A I a c t u a l ) | = | 1 1 + [ 2 R M S E ( L A I [ t ] ) σ ( Δ L A I a c t u a l ) ] 2 | ,
The uncertainty analysis for LUECanopy is complicated by the presence of high frequency components that need to be smoothed by Equation (3) in order to fully understand the temporal resolution of the retrieval. As the baseline retrieval methods with LANDSAT cannot account for the effects of the temporal smoothing because LANDSAT does not make daily measurements, the baseline retrieval must be retrained with MODIS measurements. Thus, leave-one-site-out regression is used to determine the regression coefficients in
L U E C a n o p y ¯ = p E V I 2 ¯ + q ,
where E V I 2 ¯ is the moving average of EVI2 defined as
E V I 2 ¯ [ t ] = 1 2 N 1 i = N + 1 N 1 E V I 2 [ t + i ] ,
With these leave-one-site-out regression coefficients, a baseline RMSE for the retrieval of L U E C a n o p y ¯ can be identified. In addition, as we have the benefit of a daily time series of MODIS measurements, Δ L U E C a n o p y ¯ (defined in the same way as ΔLUECanopy in Equation (2) can be determined by training a direct regression
Δ L U E C a n o p y ¯ = r ( E V I 2 ¯ [ t + Δ 1 ] E V I 2 ¯ [ t Δ + 1 ] ) + s ,
in place of using Equation (10). The regression coefficients in Equation (12) are determined by leave-one-site-out cross-validation and the performance is compared to the theoretical |r| and RMSE performance defined in Equations (8) and (9) (with LUECanopy substituted for LAI). As using Equation (12) depends on having multiple sites for cross-validation, this analysis is only performed for the actual LUECanopy measurements, while only the |r| correlation with MODIS measurements is analyzed for the modeled measurements. The analysis for LUECanopy measurements is performed between the planting and harvest dates reported for the sites; the LUECanopy analysis is not performed at US-Bi2 due to the unavailability of planting and harvest dates at this site.

2.3.5. Training LAI and LUECanopy Retrievals with HM Simulations

Lastly, in order to validate the concept of training and testing field-scale remote sensing retrievals with crop growth model simulations, we compare the performance of LAI and LUECanopy at sites other than those in Mead, Nebraska using (a) regression coefficients trained with the actual LAI and LUECanopy measurements at the Mead, Nebraska sites; and using (b) regression coefficients trained with HM modeled LAI and LUECanopy values at the Mead, Nebraska sites. These retrievals are trained and evaluated using LANDSAT measurements and the performance is reported site-by-site.

3. Results

3.1. Evaluation of HM Simulations

We first evaluate the performance of the modeled HM LAI and LUECanopy at the Mead, Nebraska sites. In Figure 1a,b, we show scatterplots between the modeled HM LAI and LUECanopy values and the actual values on the ground. As discussed in Section 2.3.1, only dates that have both LAI and LUECanopy measurements are included in Figure 1a,b for consistent comparison of the modeling performance of these two variables. The figures show strong performance for modeled LAI and LUECanopy with R2 values of 0.91 and 0.77 and RMSE values of 0.62 and 0.30, respectively; although, the bias for LUECanopy is relatively high.
In Figure 2, the performance of modeled LUECanopy and ΔLUECanopy are shown for all ground measurements of LUECanopy, not only those that also have a LAI measurement on the same date. Figure 2a shows the scatterplot of modeled LUECanopy versus actual LUECanopy with no smoothing, while Figure 2b shows the R2 value between modeled and actual LUECanopy and ΔLUECanopy at different levels of smoothing and values of Δ. As seen in Equation (3), a smoothing window of 1 represents no smoothing. Only days where modeled LUECanopy is greater than zero are included in Figure 2. In addition, a small number of days which have less than 95% of the underlying GPP time series available are not included in Figure 2.
The results in Figure 2 show that the performance of modeled LUECanopy is strong with an R2 of 0.76 in the absence of smoothing and slightly higher with smoothing. In contrast, as seen in Figure 2b, the performance of ΔLUECanopy is dependent on the level of smoothing and value of Δ, with stronger performance with longer Δ windows and more smoothing.

3.2. Regression-Based LAI and LUECanopy Retrieval

We now present the results of the retrieval of LAI and LUECanopy from LANDSAT EVI2 by Equations (4) and (5) via leave-one-site-out cross validation. In Figure 3, we present the leave-one-site-out performance for all sites combined in separate scatterplots for the LAIGROUND and FLUX datasets (prediction performed with leave-one-site-out site-by-site and then combined into a single scatter plot). Figure 3a shows the LAI retrieval scatterplot for the LAIGROUND dataset, while Figure 3b,c show the LAI and LUECanopy retrieval scatterplots for the FLUX dataset.
Figure 3 shows LAI retrieved with a R2 performance between 0.41 and 0.69 and an RMSE between 1.07 and 1.22, while LUECanopy is retrieved with an R2 performance of 0.74 and an RMSE of 0.17. In addition, the site-by-site leave-one-site-out retrieval performance and regression coefficients for the LAIGROUND dataset are shown in Table 3, while the corresponding information for the FLUX dataset is shown in Table 4. Table 3 and Table 4 also show the confidence intervals for the determined leave-one-site-out coefficients.

3.3. Satellite Retrieval and Crop Growth Model Sensitivity Analysis

We now turn to presenting the results of the crop growth model-based sensitivity analysis. First, in Figure 4, we show the RMSE of the modeled LAI with respect to the actual ground truth LAI for different simulations where three G × M parameters (the planting date, seed GDD to maturity, and planting density) are offset by various amounts from their actual values. The results in Figure 4 allow for analysis of the effect of biases in combinations of the three G × M parameters varied in the figures. The results show that with respect to the ground-truth there are several combinations of parameter bias which lead to LAI RMSEs below 0.7 against the ground-truth measurements, demonstrating ill-posedness in the inversion of LAI values to G × M parameters. As expected, the situation where none of the parameters are biased (i.e., the actual G × M parameters applied in the field, at the center of the figure), leads to a low RMSE (near 0.6), however other combinations of biases have similar RMSE. The magnitude of the error seems to be most sensitive to variations in the planting density (as seen by patterns in the variation of the performance corresponding to the frequency of the density variation); however, significant negative GDD offsets and positive planting day delays are also seen to significantly increase the error. Overall, the error is highly variable with respect to the parameter biases and many combinations of biases lead to high error (a range of LAI RMSEs from 0.6 to 1.6 is observed). This variation shows the strong sensitivity of the LAI to these three G × M inputs and the interactions between them.
In Figure 5, the sensitivity analysis from Figure 4 is reproduced with MODIS LAI retrievals instead of ground-truth LAI measurements. First, it is important to note that the analysis causes a great increase in the number of points analyzed (from N = 146 to N = 3280) and removes potential biases from a skewed distribution of growth stages as all dates are included, instead of just the dates where the ground-truth LAI measurements were taken. Secondly, the figure shows the change in modeled versus retrieved LAI error as the MODIS EVI2/LAI regression coefficients are varied. The results show the strong dependence of the error on both the regression coefficients used and the bias in the model parameters. Interestingly, although all regression coefficients show good performance for some combinations of G × M biases, some regression coefficients show significantly less sensitivity to G × M biases than others in terms of LAI error. For example, low regression slopes allow for low RMSE values at a limited number of G × M bias combinations, while high regression slopes allow for low RMSE values at a significantly greater number of G × M bias combinations. As in Figure 4, the variation in the LAI RMSE error is very sensitive to the variation of planting density, although negative GDD offsets also have a very significant effect in increasing the error. The ill-posedness of inverting the G × M factors from the MODIS measurements is seen clearly in the figure with several combinations of biases and regression coefficients leading to similar levels of LAI error. As expected, low parameter biases (near the center of the figure) lead to low LAI RMSE values, although negatively biasing the planting density appears to allow for better matchup with the MODIS measurements over a wider range of regression coefficients.

3.4. Evaluation of Uncertainty of LAI and LUECanopy Retrievals Due to Site and Growth Stage Specific Factors with Temporal Analysis

We now present the results analyzing the uncertainty of LAI and LUECanopy retrievals due to site and growth stage specific factors with temporal analysis. First, in Figure 6, we show scatterplots of retrieved versus HM modeled ΔLAI at three values of Δ (Figure 6a–c, Δ = 3, 6, 9) and compare them to the retrieval performance of HM modeled LAI itself (Figure 6d). The leave-one-out regression values from Table 3 for Mead are used to perform the retrievals. The results in Figure 6 show a rising level of performance with increasing Δ values, ranging from an R2 of 0.41 for Δ = 3 to an R2 of 0.72 at Δ = 9. The retrieval of modeled LAI itself is seen to have an R2 of 0.85 in Figure 6d.
In Table 5, we show the actual and theoretical, modeled versus retrieved |r| and RMSE for LAI itself and ΔLAI for Δ = 2 to 10. The results in Table 5 show that the actual |r| and RMSE performance of the ΔLAI retrievals significantly outperform the theoretical performance calculated with Equations (8) and (9); for example, for Δ = 4, which corresponds to a week of measurements, the actual |r| and RMSE values are 0.70 and 0.40, while the theoretical values are 0.36 and 1.46, respectively.
In Figure 7, we present the |r| correlation of the MODIS EVI2 measurements versus modeled ΔLUECanopy at different levels of smoothing and values of Δ. The results in Figure 7 show that the |r| MODIS EVI2/ΔLUECanopy correlation strongly depends on the level of smoothing and the value of Δ; however, high |r| values may be obtained when smoothing is performed at Δ values above 5.
Only days with modeled LUECanopy greater than zero are included in Figure 6 and Figure 7 and Table 5. In addition, for consistency with Figure 2, the small number of days which have less than 95% of the underlying measured GPP time series available are not included in Figure 7.
In addition to comparison of modeled values (from the Mead, Nebraska sites) in Figure 7, ΔLUECanopy retrievals are compared against the actual ΔLUECanopy measurements at all sites in the FLUX dataset. The performance against all the actual ΔLUECanopy data is shown in Figure 8 and Figure 9, while the site-by-site performance is shown in Figures S1–S10 in the Supplementary Materials.
Figure 8 and Figure 9 and Figures S1–S10 show that the actual ΔLUECanopy is retrieved with a performance at or above the theoretical performance assuming independence of retrieval error with respect to time and the relationship holds at most, but not all, sites. Poor performance at some sites, as seen in Figures S1–S10, may be explained by the large pixel size of MODIS (500 m), which can cause significant noise in the measurement of the daily time series depending on the size of the field and the inhomogeneity of the area surrounding the field [97,98]. Specifically, the mixed-pixel effect causes the signatures of neighboring pixels to be blended and it makes it difficult to separate the time-series of individual crops, especially if spring and winter crops are grown nearby [97], as is the case at some of the sites where poor performance is observed. Furthermore, the footprint of the flux tower measurements themselves depends on meteorological conditions and can be affected by process occurring outside the boundaries of the field [81,99]. Overall, however, strong performance is seen for the majority of site-years analyzed providing confidence in the retrievals.

3.5. Training LAI and LUECanopy Retrievals with HM Simulations

Lastly, we present the results indicating the performance of training the LAI and LUECanopy retrievals with HM modeled values as opposed to measured ground-truth values. In Table 6, we compare the RMSE of the LAI retrieval at sites other than Mead in the LAIGROUND dataset trained on either actual or modeled Mead LAI values, while in Table 7 we do the same for the LUECanopy retrievals in the FLUX dataset.
The results in Table 6 and Table 7 show the difference in performance in using modeled versus actual data to train the LAI regression is small, while LUECanopy retrievals perform better when trained with actual, as opposed to modeled, values.

4. Discussion

The results presented in this study outline the importance of reducing the uncertainty in the relationship between satellite measurements and crop states variables caused by site and growth stage specific factors, in particular to use using remote sensing to map the G × E × M factors affecting crop growth. The importance of reducing the uncertainty is well illustrated by Figure 5, which shows that the set of allowable G × M parameters in terms of consistency with the MODIS LAI retrievals (as measured by the RMSE) is strongly a function of the regression coefficients chosen.
Based on the “one place, one time, one equation” concept [27], the appropriate regression coefficients for each time and place are ultimately different; therefore, auxiliary information is needed to select the appropriate regression coefficient column for each site and time to retrieve G × M in Figure 5. This variability of the regression coefficients is best seen on the LAIGROUND dataset with high-resolution LANDSAT measurements in both the coefficients themselves and the large confidence intervals in Table 3, from which the range of the regression coefficients in Figure 5 was constructed. Less variability is seen on the FLUX dataset in Table 4 because this dataset has fewer points, smaller diversity in sites (points from Mead, Nebraska make up more than half the dataset), and is not designed to test the spatial variability of LAI in nearby plots in the same way as the LAIGROUND dataset; as a result, the LAIGROUND results in Table 3 are more appropriate for analyzing the variability between sites. Analyzed in conjunction with Figure 5, the regression coefficient variability in Table 3 makes it very difficult to use remote sensing for mapping G × E × M. This is because, as illustrated in Figure 4 and Figure 5, the retrieval of G × E × M is difficult due to equifinality (i.e., “multiple combinations of parameters leading to similar simulation accuracy”) [13] and, especially when the observations are uncertain remote sensing retrievals, is ill-posed. Figure 4 does a good job of showing the ill-posedness of the G × M retrieval even when using ground-truth LAI measurements; interestingly, due to the availability of the entire time series when using MODIS measurements in Figure 5, some combinations of G × M identified as probable in Figure 4 are not probable in Figure 5 for any combination of regression coefficients. This is an illustration of the importance of the number of measurements [13,100] needed to perform G × E × M retrievals and the frequent, low-cost observations provided by satellites may be one of the most promising technologies to achieve that goal [101].
Although the uncertainty caused by site and growth stage specific secondary factors is well-known [27,28,29,50,51], it is difficult to isolate it from other sources. One approach to understand it [29,50,51] is to include variables connected to the secondary factors that cause it in the regression methodology. Unfortunately, this approach requires that the secondary factors causing the uncertainty are known and recorded or measured prior to the analysis being conducted. As a result, these studies can miss some of the factors causing the issue and underestimate its extent. Another approach is to train a global relationship between the satellite variables and crop state variables, ignoring the secondary factors [27]. In this case [27], the issue is seen from the variability of the regression coefficients, as in our analysis in Table 3 and Table 4, as well as indirectly from the variability in the leave-one-site-out RMSE error. However, this method cannot be used to exclude other sources of uncertainty from the retrieval, such as random error and the mixed-pixel effect [102].
In this study, temporal analysis is used to avoid these alternate sources [27,29,50,51] of uncertainty in determining the portion caused by site and growth stage specific secondary factors. The results in Table 5 show that the modeled ΔLAI is retrieved from the MODIS measurements with significantly better performance in terms of both |r| and RMSE as compared to the theoretical values assuming temporal independence of error, indicating significant site and time correlation of error. These results are also reproduced with actual ΔLUECanopy measurements across multiple CO2 flux tower sites in Figure 8 and Figure 9. Both the results with modeled ΔLAI and actual ΔLUECanopy indicate a significant portion of the error can be removed by either predicting the secondary factors [29,50,51] or developing better methods to remove their influence, such as identifying vegetation indices less sensitive to the secondary factors [32,45,46]. The difference between the actual and theoretical |r| and RMSE for both ΔLAI and ΔLUECanopy provide an indication of the possible reduction in uncertainty by addressing the issue with secondary factors. The change in performance with respect to the value of Δ is driven by two factors:
  • As Δ increases, the correlation between the error in the retrieved LAI or LUECanopy at t2 relative to t1 decreases because the measurements are more likely to be in different growth stages.
  • As Δ increases, the magnitude of the retrieved ΔLAI or ΔLUECanopy increases relative to the remaining error which is not cancelled when calculating the change in the retrieved variables from the variables themselves, i.e., e[t2] − e[t1].
As a result of these opposing error-influencing forces, a single value for the improvement in the performance that could be obtained by reducing the influence of the secondary factors cannot be reported; however, as seen from Table 5 and Figure 8 and Figure 9 the improvement can be quite dramatic. For example, for ΔLAI, the actual |r| at Δ = 2 is 0.52 (compared to a theoretical value of 0.13), while the actual |r| at Δ = 5 is 0.75 (compared to a theoretical value of 0.45).
Furthermore, the retrieval of ΔLAI and ΔLUECanopy is also useful as a measure of the timescale of the sensitivity of the MODIS measurements to changes in the canopy structure and crop status. Good responsiveness to time-sensitive processes is important for several applications of crop remote sensing. For example, good responsiveness is important in monitoring phenology/crop growth stage [103,104,105], in-season detection of nitrogen [106,107], water [107], and disease [108] stresses, and measurement of change in canopy structure during important growth stages [109]; these applications have proven useful in crop growth modeling [59], precision agriculture [110], and phenotyping for breeding selection [109], respectively. Our results show that satellite measurements can be used to detect changes in LAI and LUECanopy faster and with higher accuracy than would be expected if the error in LAI and LUECanopy retrievals were not autocorrelated in time. As a result, we also show the potential to rapidly detect growth and stress related changes in crop state variables with greater precision than that would be inferred from looking at generic performance validation studies [27,28].
The analysis used in this study relied on strong crop growth model simulation performance to expand the dataset of ground-truth LAI values to daily resolution. The strong performance of the HM simulations at the Mead, Nebraska sites, seen in [90] and Figure 1 and Figure 2, provides a potential path [67] for future research to expand the development of testing agronomic satellite retrievals to a wide variety of G x E x M factors with farmer-provided agromanagement data. The results in Table 6 show that using HM simulation data from Mead, Nebraska to train LAI retrievals can provide nearly identical performance to using actual ground-truth LAI measurements from Mead, while Table 7 shows there are some relatively significant biases in using modeled LUECanopy to perform the training. The results for training LAI retrievals on HM simulation data show the potential of using farmer-provided agromanagement data to train, test, and improve retrieval algorithms, although a significantly greater number of sites is needed to understand the generalizability and biases in this approach. Nevertheless, the potential of reducing the uncertainty in the retrieval of crop state variables and the potential to map G × E × M factors shown in this study provides strong support for pursuing this collocated crop growth model simulation approach in agricultural remote sensing and should encourage researchers to increase their collaborative efforts with farmers [68,111].

5. Conclusions

Overall, this study showed that the uncertainty in the relationship between satellite measurements and crop state variables caused by site and growth stage specific factors is significant and that this uncertainty leads to significant difficulties in using remotely sensed data to retrieve the genotype × environment × management (G × E × M) factors affecting crop growth. Specifically, we performed an extensive temporal analysis and retrieved the temporal change in the state variables to show the amount of uncertainty caused by this secondary factor variability. We also conducted a joint sensitivity analysis of the remote sensing regression parameters and crop model genotype x management (G × M) parameters to illustrate the ill-posedness of retrieving G × E × M factors from satellite measurements. This analysis demonstrated the criticalness of reducing the uncertainty in the relationship between satellite measurements and crop state variables to make the retrieval more feasible. The study shows the need for further data collection and model development that can ultimately lead to methods that will minimize the secondary uncertainty caused by site and growth stage specific factors. In addition, further work needs to be conducted to address the application of the methods to use training data in regions where biotic stresses are poorly controlled and where, unlike the case in highly developed commercial agriculture systems [7], crop growth models show significant uncertainties in predicting actual yields (as opposed to potential yields) due to suboptimal management [112]. This research is critical to achieving the goal of mapping G × E × M factors on a global scale, which can improve our ability to make predictions about the global agricultural system [113].

Supplementary Materials

The following are available online at https://www.mdpi.com/2072-4292/11/16/1928/s1, Figures S1–S10: Actual versus theoretical |r| and RMSE performance of the ΔLUECanopy retrievals site-by-site at different levels of smoothing and values of Δ.

Author Contributions

Conceptualization, N.L., Y.K., P.C., F.M., and B.G.; Methodology, N.L. and B.G.; Software, N.L.; Validation, N.L.; Formal analysis, N.L.; Investigation, N.L., Y.K., M.O., and V.M.; Resources, Y.K., M.O., and V.M.; Data curation, N.L., Y.K., M.O., and V.M.; Writing—original draft preparation, N.L.; Writing—review and editing, N.L., Y.K., M.O., V.M., P.C., F.M., and B.G.; Visualization, N.L.; Supervision, F.M. and B.G.; Project administration, F.M. and B.G.

Funding

This work was supported by NASA Headquarters under the NASA Earth and Space Science Fellowship Program—Grant 80NSSC17K0339. We would also like to acknowledge partial support by NOAA under the grant CREST agreement #NA16SEC4810008. In addition, data used in this study was funded by a variety of sources. Observations collected at the US Southern Great Plains site (US-ARM) were supported by the Office of Biological and Ecological Research of the US Department of Energy under contract no. DE-AC02-05ch11231 as part of the Atmospheric Radiation Measurement (ARM) Program. Observations at DE-Kli were funded by EU project CARBOEUROPE-IP, by German BMBF project ICOS-D and by the state of Saxony (TU Dresden, LfULG). Data acquisition at FR-Lam are mainly funded by the Institut National des Sciences de l’Univers (INSU) through the ICOS and OSR SW observatories. Facilities and staff at FR-Lam are funded and supported by the Observatory Midi-Pyrenean, the University Paul Sabatier of Toulouse, CNRS (Centre National de la Recherche Scientifique), CNES (Centre National d’Etude Spatial), and IRD (Institut de Recherche pour le Développement). Facilities and staff at FR-Gri were funded and supported by INRA, AgroParisTech, European projects NitroEurope-IP (EU-) ECLAIRE-FP7 (EU-), INGOS-FP7, Accent, and ESF-NinE, French National Research Agency (ANR) projects Vulnoz, CNRS, Ile de France region R2DS, CNRS-INSU, Fr-De Photona, Programme LabEx BASC, ANR-11-LABX-0034.

Acknowledgments

We are grateful to all of the investigators and personnel who collected and provided the ground truth data used in this study. In particularly, we would to like to acknowledge the following investigators and personnel (site identification in parenthesis): Andy Suyker (US-Ne1, US-Ne2, and US-Ne3), John Baker (US-Ro1), Tim Griffis (US-Ro1), Cody Winker (US-Ro1), Dennis Baldocchi (US-Bi2), Daphne Szutu (US-Bi2), Sebastien C. Biraud (US-ARM), Christian Bernhofer (DE-Kli), Thomas Grünwald (DE-Kli), Benjamin Loubet (FR-Gri), Pierre Cellier (FR-Gri), Pauline Buysse (FR-Gri), Aurore Brut (FR-Lam—data processing), Tiphaine Tallec (FR-Lam—data processing), Franck Granouillac (FR-Lam—technical support in field and site management), Bartosz Zawilski (FR-Lam—technical support in field and site management), Vincenzo Magliulo (IT-BCi), Bart Kruijt (NL-Lan), and Wilma Jans (NL-Lan). In addition, we want to thank Mr. Andréoni, farmer, and the Ecole Supérieure d’Agriculture de Purpan for accommodating the measurement devices in their field at the FR-Lam site. Furthermore, we want to thank Dominique Tristant, farmer, at the AgroParisTech Farm for lending his field to the FR-Gri station. This work is part of a dissertation to be submitted to The City College of New York entitled “Improving Retrievals of Crop Vegetation Parameters from Remote Sensing Data”.

Conflicts of Interest

The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

References

  1. Teixeira, E.I.; Zhao, G.; de Ruiter, J.; Brown, H.; Ausseil, A.-G.; Meenken, E.; Ewert, F. The interactions between genotype, management and environment in regional crop modeling. Eur. J. Agron. 2017, 88, 106–115. [Google Scholar] [CrossRef]
  2. Archontoulis, S.V.; Miguez, F.E.; Moore, K.J. Evaluating APSIM maize, soil water, soil nitrogen, manure, and soil temperature modules in the Midwestern United States. Agron. J. 2014, 106, 1025–1040. [Google Scholar] [CrossRef]
  3. Keating, B.; Carberry, P.; Hammer, G.; Probert, M.; Robertson, M.; Holzworth, D.; Huth, N.; Hargreaves, J.N.; Meinke, H.; Hochman, Z.; et al. An overview of APSIM, a model designed for farming systems simulation. Eur. J. Agron. 2003, 18, 267–288. [Google Scholar] [CrossRef] [Green Version]
  4. Donatelli, M.; Magarey, R.D.; Bregaglio, S.; Willocquet, L.; Whish, J.P.M.; Savary, S. Modeling the impacts of pests and diseases on agricultural systems. Agric. Syst. 2017, 155, 213–224. [Google Scholar] [CrossRef]
  5. Grassini, P.; Thorburn, J.; Burr, C.; Cassman, K.G. High-yield irrigated maize in the Western, U.S. Corn Belt: I. On-farm yield, yield potential, and impact of agronomic practices. Field Crops Res. 2011, 120, 142–150. [Google Scholar] [CrossRef]
  6. Morell, F.J.; Yang, H.S.; Cassman, K.G.; Wart, J.V.; Elmore, R.W.; Licht, M.; Coulter, J.A.; Ciampitti, I.A.; Pittelkow, C.M.; Brouder, S.M.; et al. Can crop simulation models be used to predict local to regional maize yields and total production in the U.S. Corn Belt? Field Crops Res. 2016, 192, 1–12. [Google Scholar] [CrossRef] [Green Version]
  7. Carberry, P.S.; Hochman, Z.; Hunt, J.R.; Dalgliesh, N.P.; McCown, R.L.; Whish, J.P.M.; Robertson, M.J.; Foale, M.A.; Poulton, P.L.; van Rees, H. Re-inventing model-based decision support with Australian dryland farmers. 3. Relevance of APSIM to commercial crops. Crop Pasture Sci. 2009, 60, 1044–1056. [Google Scholar] [CrossRef] [Green Version]
  8. Jin, X.; Kumar, L.; Li, Z.; Feng, H.; Xu, X.; Yang, G.; Wang, J. A review of data assimilation of remote sensing and crop models. Eur. J. Agron. 2018, 92, 141–152. [Google Scholar] [CrossRef]
  9. Dorigo, W.A.; Zurita-Milla, R.; de Wit, A.J.W.; Brazile, J.; Singh, R.; Schaepman, M.E. A review on reflective remote sensing and data assimilation techniques for enhanced agroecosystem modeling. Int. J. Appl. Earth Obs. Geoinf. 2007, 9, 165–193. [Google Scholar] [CrossRef]
  10. Guérif, M.; Duke, C. Adjustment procedures of a crop model to the site specific characteristics of soil and crop using remote sensing data assimilation. Agric. Ecosyst. Environ. 2000, 81, 57–69. [Google Scholar] [CrossRef]
  11. Constantin, J.; Raynal, H.; Casellas, E.; Hoffmann, H.; Bindi, M.; Doro, L.; Eckersten, H.; Gaiser, T.; Grosz, B.; Haas, E.; et al. Management and spatial resolution effects on yield and water balance at regional scale in crop models. Agric. For. Meteorol. 2019, 275, 184–195. [Google Scholar] [CrossRef]
  12. He, J.; Jones, J.W.; Graham, W.D.; Dukes, M.D. Influence of likelihood function choice for estimating crop model parameters using the generalized likelihood uncertainty estimation method. Agric. Syst. 2010, 103, 256–264. [Google Scholar] [CrossRef]
  13. He, D.; Wang, E.; Wang, J.; Robertson, M.J. Data requirement for effective calibration of process-based crop models. Agric. For. Meteorol. 2017, 234–235, 136–148. [Google Scholar] [CrossRef]
  14. Confalonieri, R.; Bregaglio, S.; Acutis, M. Quantifying uncertainty in crop model predictions due to the uncertainty in the observations used for calibration. Ecol. Model. 2016, 328, 72–77. [Google Scholar] [CrossRef]
  15. Kang, Y.; Özdoğan, M. Field-level crop yield mapping with Landsat using a hierarchical data assimilation approach. Remote Sens. Environ. 2019, 228, 144–163. [Google Scholar] [CrossRef]
  16. Battude, M.; Al Bitar, A.; Morin, D.; Cros, J.; Huc, M.; Marais Sicre, C.; Le Dantec, V.; Demarez, V. Estimating maize biomass and yield over large areas using high spatial and temporal resolution Sentinel-2 like remote sensing data. Remote Sens. Environ. 2016, 184, 668–681. [Google Scholar] [CrossRef]
  17. Ines, A.V.M.; Das, N.N.; Hansen, J.W.; Njoku, E.G. Assimilation of remotely sensed soil moisture and vegetation with a crop simulation model for maize yield prediction. Remote Sens. Environ. 2013, 138, 149–164. [Google Scholar] [CrossRef] [Green Version]
  18. Jégo, G.; Pattey, E.; Liu, J. Using Leaf Area Index, retrieved from optical imagery, in the STICS crop model for predicting yield and biomass of field crops. Field Crops Res. 2012, 131, 63–74. [Google Scholar] [CrossRef]
  19. Battisti, R.; Sentelhas, P.C.; Boote, K.J. Inter-comparison of performance of soybean crop simulation models and their ensemble in southern Brazil. Field Crops Res. 2017, 200, 28–37. [Google Scholar] [CrossRef]
  20. Angulo, C.; Rötter, R.; Lock, R.; Enders, A.; Fronzek, S.; Ewert, F. Implication of crop model calibration strategies for assessing regional impacts of climate change in Europe. Agric. For. Meteorol. 2013, 170, 32–46. [Google Scholar] [CrossRef]
  21. Rosenzweig, C.; Elliott, J.; Deryng, D.; Ruane, A.C.; Müller, C.; Arneth, A.; Boote, K.J.; Folberth, C.; Glotter, M.; Khabarov, N.; et al. Assessing agricultural risks of climate change in the 21st century in a global gridded crop model intercomparison. Proc. Natl. Acad. Sci. USA 2014, 111, 3268–3273. [Google Scholar] [CrossRef]
  22. Lecerf, R.; Ceglar, A.; López-Lozano, R.; Van Der Velde, M.; Baruth, B. Assessing the information in crop model and meteorological indicators to forecast crop yield over Europe. Agric. Syst. 2019, 168, 191–202. [Google Scholar] [CrossRef]
  23. Woodard, J.D. Integrating high resolution soil data into federal crop insurance policy: Implications for policy and conservation. Environ. Sci. Policy 2016, 66, 93–100. [Google Scholar] [CrossRef] [Green Version]
  24. Grassini, P.; Wolf, J.; Tittonell, P.; Hochman, Z. Yield gap analysis with local to global relevance—A review. Field Crop. Res. 2013, 143, 4–17. [Google Scholar] [CrossRef]
  25. Ramirez-Villegas, J.; Koehler, A.-K.; Challinor, A.J. Assessing uncertainty and complexity in regional-scale crop model simulations. Eur. J. Agron. 2017, 88, 84–95. [Google Scholar] [CrossRef] [Green Version]
  26. Reynolds, M.; Kropff, M.; Crossa, J.; Koo, J.; Kruseman, G.; Molero, A.M.; Rutkoski, J.; Schulthess, U.; Balwinder-Singh; Sonder, K.; et al. Role of Modeling in International Crop Research: Overview and Some Case Studies. Agronomy 2018, 8, 291. [Google Scholar] [CrossRef]
  27. Kang, Y.; Özdoğan, M.; Zipper, S.; Román, M.; Walker, J.; Hong, S.; Marshall, M.; Magliulo, V.; Moreno, J.; Alonso, L.; et al. How Universal Is the Relationship between Remotely Sensed Vegetation Indices and Crop Leaf Area Index? A Global Assessment. Remote Sens. 2016, 8, 597. [Google Scholar] [CrossRef]
  28. Corti, M.; Cavalli, D.; Cabassi, G.; Marino Gallina, P.; Bechini, L. Does remote and proximal optical sensing successfully estimate maize variables? A review. Eur. J. Agron. 2018, 99, 37–50. [Google Scholar] [CrossRef]
  29. Wang, Y.; Zhang, K.; Tang, C.; Cao, Q.; Tian, Y.; Zhu, Y.; Cao, W.; Liu, X.; Wang, Y.; Zhang, K.; et al. Estimation of Rice Growth Parameters Based on Linear Mixed-Effect Model Using Multispectral Images from Fixed-Wing Unmanned Aerial Vehicles. Remote Sens. 2019, 11, 1371. [Google Scholar] [CrossRef]
  30. Clevers, J.; Kooistra, L.; van den Brande, M.; Clevers, J.G.P.W.; Kooistra, L.; Van den Brande, M.M.M. Using Sentinel-2 Data for Retrieving LAI and Leaf and Canopy Chlorophyll Content of a Potato Crop. Remote Sens. 2017, 9, 405. [Google Scholar] [CrossRef]
  31. Gonsamo, A.; Chen, J.M. Improved LAI Algorithm Implementation to MODIS Data by Incorporating Background, Topography, and Foliage Clumping Information. IEEE Trans. Geosci. Remote Sens. 2014, 52, 1076–1088. [Google Scholar] [CrossRef]
  32. Liu, J.; Pattey, E.; Jégo, G. Assessment of vegetation indices for regional crop green LAI estimation from Landsat images over multiple growing seasons. Remote Sens. Environ. 2012, 123, 347–358. [Google Scholar] [CrossRef]
  33. Li, Z.; Jin, X.; Yang, G.; Drummond, J.; Yang, H.; Clark, B.; Li, Z.; Zhao, C.; Li, Z.; Jin, X.; et al. Remote Sensing of Leaf and Canopy Nitrogen Status in Winter Wheat (Triticum aestivum L.) Based on N-PROSAIL Model. Remote Sens. 2018, 10, 1463. [Google Scholar] [CrossRef]
  34. Boren, E.J.; Boschetti, L.; Johnson, D.M.; Boren, E.J.; Boschetti, L.; Johnson, D.M. Characterizing the Variability of the Structure Parameter in the PROSPECT Leaf Optical Properties Model. Remote Sens. 2019, 11, 1236. [Google Scholar] [CrossRef]
  35. Gitelson, A.A.; Peng, Y.; Arkebauer, T.J.; Schepers, J. Relationships between gross primary production, green LAI, and canopy chlorophyll content in maize: Implications for remote sensing of primary production. Remote Sens. Environ. 2014, 144, 65–72. [Google Scholar] [CrossRef]
  36. Pinter, P.J.; Jackson, R.D.; Ezra, C.E.; Gausman, H.W. Sun-angle and canopy-architecture effects on the spectral reflectance of six wheat cultivars. Int. J. Remote Sens. 1985, 6, 1813–1825. [Google Scholar] [CrossRef]
  37. Jacquemoud, S.; Verhoef, W.; Baret, F.; Bacour, C.; Zarco-Tejada, P.J.; Asner, G.P.; François, C. PROSPECT + SAIL models: A review of use for vegetation characterization. Remote Sens. Environ. 2009, 113, S56–S66. [Google Scholar] [CrossRef]
  38. Combal, B.; Baret, F.; Weiss, M.; Trubuil, A.; Macé, D.; Pragnère, A.; Myneni, R.; Knyazikhin, Y.; Wang, L. Retrieval of canopy biophysical variables from bidirectional reflectance: Using prior information to solve the ill-posed inverse problem. Remote Sens. Environ. 2003, 84, 1–15. [Google Scholar] [CrossRef]
  39. Houborg, R.; McCabe, M.F. Adapting a regularized canopy reflectance model (REGFLEC) for the retrieval challenges of dryland agricultural systems. Remote Sens. Environ. 2016, 186, 105–120. [Google Scholar] [CrossRef] [Green Version]
  40. Houborg, R.; McCabe, M.; Cescatti, A.; Gao, F.; Schull, M.; Gitelson, A. Joint leaf chlorophyll content and leaf area index retrieval from Landsat data using a regularized model inversion system (REGFLEC). Remote Sens. Environ. 2015, 159, 203–221. [Google Scholar] [CrossRef] [Green Version]
  41. Xiao, Z.; Liang, S.; Wang, J.; Song, J.; Wu, X. A Temporally Integrated Inversion Method for Estimating Leaf Area Index From MODIS Data. IEEE Trans. Geosci. Remote Sens. 2009, 47, 2536–2545. [Google Scholar] [CrossRef]
  42. Koetz, B.; Baret, F.; Poilvé, H.; Hill, J. Use of coupled canopy structure dynamic and radiative transfer models to estimate biophysical canopy characteristics. Remote Sens. Environ. 2005, 95, 115–124. [Google Scholar] [CrossRef]
  43. Atzberger, C. Object-based retrieval of biophysical canopy variables using artificial neural nets and radiative transfer models. Remote Sens. Environ. 2004, 93, 53–67. [Google Scholar] [CrossRef]
  44. Dorigo, W.; Richter, R.; Baret, F.; Bamler, R.; Wagner, W.; Dorigo, W.; Richter, R.; Baret, F.; Bamler, R.; Wagner, W. Enhanced automated canopy characterization from hyperspectral data by a novel two step radiative transfer model inversion approach. Remote Sens. 2009, 1, 1139–1170. [Google Scholar] [CrossRef]
  45. Jin, X.; Li, Z.; Feng, H.; Xu, X.; Yang, G. Newly Combined Spectral Indices to Improve Estimation of Total Leaf Chlorophyll Content in Cotton. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2014, 7, 4589–4600. [Google Scholar] [CrossRef]
  46. Xiao, Y.; Zhao, W.; Zhou, D.; Gong, H. Sensitivity Analysis of Vegetation Reflectance to Biochemical and Biophysical Variables at Leaf, Canopy, and Regional Scales. IEEE Trans. Geosci. Remote Sens. 2014, 52, 4014–4024. [Google Scholar] [CrossRef]
  47. Verrelst, J.; Camps-Valls, G.; Muñoz-Marí, J.; Rivera, J.P.; Veroustraete, F.; Clevers, J.G.P.W.; Moreno, J. Optical remote sensing and the retrieval of terrestrial vegetation bio-geophysical properties—A review. ISPRS J. Photogramm. Remote Sens. 2015, 108, 273–290. [Google Scholar] [CrossRef]
  48. Darvishzadeh, R.; Atzberger, C.; Skidmore, A.; Schlerf, M. Mapping grassland leaf area index with airborne hyperspectral imagery: A comparison study of statistical approaches and inversion of radiative transfer models. ISPRS J. Photogramm. Remote Sens. 2011, 66, 894–906. [Google Scholar] [CrossRef]
  49. Atzberger, C.; Guérif, M.; Baret, F.; Werner, W. Comparative analysis of three chemometric techniques for the spectroradiometric assessment of canopy chlorophyll content in winter wheat. Comput. Electron. Agric. 2010, 73, 165–173. [Google Scholar] [CrossRef]
  50. Kiala, Z.; Odindi, J.; Mutanga, O.; Peerbhay, K. Comparison of partial least squares and support vector regressions for predicting leaf area index on a tropical grassland using hyperspectral data. J. Appl. Remote Sens. 2016, 10, 036015. [Google Scholar] [CrossRef]
  51. Wang, L.; Chang, Q.; Li, F.; Yan, L.; Huang, Y.; Wang, Q.; Luo, L.; Wang, L.; Chang, Q.; Li, F.; et al. Effects of growth stage development on paddy rice leaf area index prediction models. Remote Sens. 2019, 11, 361. [Google Scholar] [CrossRef]
  52. Machwitz, M.; Giustarini, L.; Bossung, C.; Frantz, D.; Schlerf, M.; Lilienthal, H.; Wandera, L.; Matgen, P.; Hoffmann, L.; Udelhoven, T. Enhanced biomass prediction by assimilating satellite data into a crop growth model. Environ. Model. Softw. 2014, 62, 437–453. [Google Scholar] [CrossRef]
  53. Weiss, M.; Troufleau, D.; Baret, F.; Chauki, H.; Prévot, L.; Olioso, A.; Bruguier, N.; Brisson, N. Coupling canopy functioning and radiative transfer models for remote sensing data assimilation. Agric. For. Meteorol. 2001, 108, 113–128. [Google Scholar] [CrossRef]
  54. Zhang, L.; Guo, C.L.; Zhao, L.Y.; Zhu, Y.; Cao, W.X.; Tian, Y.C.; Cheng, T.; Wang, X. Estimating wheat yield by integrating the WheatGrow and PROSAIL models. Field Crops Res. 2016, 192, 55–66. [Google Scholar] [CrossRef]
  55. Thorp, K.R.; Wang, G.; West, A.L.; Moran, M.S.; Bronson, K.F.; White, J.W.; Mon, J. Estimating crop biophysical properties from remote sensing data by inverting linked radiative transfer and ecophysiological models. Remote Sens. Environ. 2012, 124, 224–233. [Google Scholar] [CrossRef] [Green Version]
  56. Jin, Z.; Prasad, R.; Shriver, J.; Zhuang, Q. Crop model- and satellite imagery-based recommendation tool for variable rate N fertilizer application for the US Corn system. Precis. Agric. 2017, 18, 779–800. [Google Scholar] [CrossRef]
  57. Xiong, W.; Skalský, R.; Porter, C.H.; Balkovič, J.; Jones, J.W.; Yang, D. Calibration-induced uncertainty of the EPIC model to estimate climate change impact on global maize yield. J. Adv. Model. Earth Syst. 2016, 8, 1358–1375. [Google Scholar] [CrossRef] [Green Version]
  58. Tatsumi, K. Effects of automatic multi-objective optimization of crop models on corn yield reproducibility in the U.S.A. Ecol. Model. 2016, 322, 124–137. [Google Scholar] [CrossRef]
  59. Jägermeyr, J.; Frieler, K. Spatial variations in crop growing seasons pivotal to reproduce global fluctuations in maize and wheat yields. Sci. Adv. 2018, 4, eaat4517. [Google Scholar] [CrossRef] [Green Version]
  60. Müller, C.; Elliott, J.; Chryssanthacopoulos, J.; Arneth, A.; Balkovic, J.; Ciais, P.; Deryng, D.; Folberth, C.; Glotter, M.; Hoek, S.; et al. Global gridded crop model evaluation: Benchmarking, skills, deficiencies and implications. Geosci. Model Dev. 2017, 10, 1403–1422. [Google Scholar] [CrossRef]
  61. Challinor, A.; Martre, P.; Asseng, S.; Thornton, P.; Ewert, F. Making the most of climate impacts ensembles. Nat. Clim. Chang. 2014, 4, 77–80. [Google Scholar] [CrossRef] [Green Version]
  62. Soufizadeh, S.; Munaro, E.; McLean, G.; Massignam, A.; van Oosterom, E.J.; Chapman, S.C.; Messina, C.; Cooper, M.; Hammer, G.L. Modeling the nitrogen dynamics of maize crops—Enhancing the APSIM maize model. Eur. J. Agron. 2018, 100, 118–131. [Google Scholar] [CrossRef]
  63. Heng, L.K.; Hsiao, T.; Evett, S.; Howell, T.; Steduto, P. Validating the FAO AquaCrop Model for Irrigated and Water Deficient Field Maize. Agron. J. 2009, 101, 488. [Google Scholar] [CrossRef]
  64. Yang, H.; Dobermann, A.; Lindquist, J.; Walters, D.; Arkebauer, T.; Cassman, K. Hybrid-Maize—A maize simulation model that combines two crop modeling approaches. Field Crops Res. 2004, 87, 131–154. [Google Scholar] [CrossRef]
  65. Ben Nouna, B.; Katerji, N.; Mastrorilli, M. Using the CERES-Maize model in a semi-arid Mediterranean environment. Evaluation of model performance. Eur. J. Agron. 2000, 13, 309–322. [Google Scholar] [CrossRef]
  66. Aggarwal, P.K. Uncertainties in crop, soil and weather inputs used in growth models: Implications for simulated outputs and their applications. Agric. Syst. 1995, 48, 361–384. [Google Scholar] [CrossRef]
  67. Levitan, N.; Gross, B. Utilizing Collocated Crop Growth Model Simulations to Train Agronomic Satellite Retrieval Algorithms. Remote Sens. 2018, 10, 1968. [Google Scholar] [CrossRef]
  68. Wolfert, S.; Ge, L.; Verdouw, C.; Bogaardt, M.-J. Big Data in Smart Farming—A review. Agric. Syst. 2017, 153, 69–80. [Google Scholar] [CrossRef]
  69. Baldocchi, D.; Falge, E.; Gu, L.; Olson, R.; Hollinger, D.; Running, S.; Anthoni, P.; Bernhofer, C.; Davis, K.; Evans, R.; et al. FLUXNET: A New Tool to Study the Temporal and Spatial Variability of Ecosystem–Scale Carbon Dioxide, Water Vapor, and Energy Flux Densities. Bull. Am. Meteorol. Soc. 2001, 82, 2415–2434. [Google Scholar] [CrossRef]
  70. Pastorello, G.; Papale, D.; Chu, H.; Trotta, C.; Agarwal, D.; Canfora, E.; Baldocchi, D.; Torn, M. A new data set to keep a sharper eye on land-air exchanges. Eos, 17 April 2017. [Google Scholar] [CrossRef]
  71. Reichstein, M.; Falge, E.; Baldocchi, D.; Papale, D.; Aubinet, M.; Berbigier, P.; Bernhofer, C.; Buchmann, N.; Gilmanov, T.; Granier, A.; et al. On the separation of net ecosystem exchange into assimilation and ecosystem respiration: Review and improved algorithm. Glob. Chang. Biol. 2005, 11, 1424–1439. [Google Scholar] [CrossRef]
  72. Wutzler, T.; Lucas-Moffat, A.; Migliavacca, M.; Knauer, J.; Sickel, K.; Šigut, L.; Menzer, O.; Reichstein, M. Basic and extensible post-processing of eddy covariance flux data with REddyProc. Biogeosciences 2018, 15, 5015–5030. [Google Scholar] [CrossRef] [Green Version]
  73. Jonckheere, I.; Fleck, S.; Nackaerts, K.; Muys, B.; Coppin, P.; Weiss, M.; Baret, F. Review of methods for in situ leaf area index determination: Part, I. Theories, sensors and hemispherical photography. Agric. For. Meteorol. 2004, 121, 19–35. [Google Scholar] [CrossRef]
  74. Schmidt, G.; Jenkerson, C.B.; Masek, J.; Vermote, E.; Gao, F. Landsat Ecosystem Disturbance Adaptive Processing System (LEDAPS) Algorithm Description; U.S. Geological Survey: Reston, VA, USA, 2013.
  75. Schaaf, C.; Wang, Z. MODIS/Terra and Aqua Nadir BRDF-Adjusted Reflectance Daily L3 Global 500 m SIN Grid V006; NASA EOSDIS Land Processes DAAC: Sioux Falls, SD, USA, 2015.
  76. Hargrove, W.W.; Hoffman, F.M.; Law, B.E. New analysis reveals representativeness of the AmeriFlux network. EOS Trans. Am. Geophys. Union 2003, 84, 529–535. [Google Scholar] [CrossRef]
  77. Griffis, T.J.; Baker, J.M.; Zhang, J. Seasonal dynamics and partitioning of isotopic CO2 exchange in a C3/C4 managed ecosystem. Agric. For. Meteorol. 2005, 132, 1–19. [Google Scholar] [CrossRef]
  78. Hemes, K.S.; Chamberlain, S.D.; Eichelmann, E.; Anthony, T.; Valach, A.; Kasak, K.; Szutu, D.; Verfaillie, J.; Silver, W.L.; Baldocchi, D.D. Assessing the carbon and climate benefit of restoring degraded agricultural peat soils to managed wetlands. Agric. For. Meteorol. 2019, 268, 202–214. [Google Scholar] [CrossRef]
  79. Raz-Yaseef, N.; Billesbach, D.P.; Fischer, M.L.; Biraud, S.C.; Gunter, S.A.; Bradford, J.A.; Torn, M.S. Vulnerability of crops and native grasses to summer drying in the U.S. Southern Great Plains. Agric. Ecosyst. Environ. 2015, 213, 209–218. [Google Scholar] [CrossRef] [Green Version]
  80. Prescher, A.-K.; Grünwald, T.; Bernhofer, C. Land use regulates carbon budgets in eastern Germany: From NEE to NBP. Agric. For. Meteorol. 2010, 150, 1016–1025. [Google Scholar] [CrossRef]
  81. Loubet, B.; Laville, P.; Lehuger, S.; Larmanou, E.; Fléchard, C.; Mascher, N.; Genermont, S.; Roche, R.; Ferrara, R.M.; Stella, P.; et al. Carbon, nitrogen and Greenhouse gases budgets over a four years crop rotation in northern France. Plant Soil 2011, 343, 109–137. [Google Scholar] [CrossRef]
  82. Béziat, P.; Ceschia, E.; Dedieu, G. Carbon balance of a three crop succession over two cropland sites in South West France. Agric. For. Meteorol. 2009, 149, 1628–1645. [Google Scholar] [CrossRef] [Green Version]
  83. Vitale, L.; Di Tommasi, P.; Arena, C.; Fierro, A.; Virzo De Santo, A.; Magliulo, V. Effects of water stress on gas exchange of field grown Zea mays L. in Southern Italy: An analysis at canopy and leaf level. Acta Physiol. Plant. 2007, 29, 317–326. [Google Scholar] [CrossRef]
  84. Moors, E.J.; Jacobs, C.; Jans, W.; Supit, I.; Kutsch, W.L.; Bernhofer, C.; Béziat, P.; Buchmann, N.; Carrara, A.; Ceschia, E.; et al. Variability in carbon exchange of European croplands. Agric. Ecosyst. Environ. 2010, 139, 325–335. [Google Scholar] [CrossRef]
  85. European Space Agency Earth Observation Campaigns Data. Available online: https://earth.esa.int/web/guest/campaigns (accessed on 2 September 2013).
  86. Marshall, M.; Thenkabail, P. Developing in situ non-destructive estimates of crop biomass to address issues of scale in remote sensing. Remote Sens. 2015, 7, 808–835. [Google Scholar] [CrossRef]
  87. Kim, S.-H.; Hong, S.Y.; Sudduth, K.A.; Kim, Y.; Lee, K. Comparing LAI Estimates of Corn and Soybean from Vegetation Indices of Multi-resolution Satellite Images. Korean J. Remote Sens. 2012, 28, 597–609. [Google Scholar] [CrossRef] [Green Version]
  88. Merlin, O.; Walker, J.P.; Kalma, J.D.; Kim, E.J.; Hacker, J.; Panciera, R.; Young, R.; Summerell, G.; Hornbuckle, J.; Hafeez, M.; et al. The NAFE’06 data set: Towards soil moisture retrieval at intermediate resolution. Adv. Water Resour. 2008, 31, 1444–1455. [Google Scholar] [CrossRef]
  89. Anderson, M. SMEX02 Regional Vegetation Sampling Data, Iowa; National Snow and Ice Data Center: Boulder, CO, USA, 2003. [Google Scholar]
  90. Zhan, M.; Liska, A.J.; Nguy-Robertson, A.L.; Suyker, A.E.; Pelton, M.P.; Yang, H. Modeled and Measured Ecosystem Respiration in Maize–Soybean Systems Over 10 Years. Agron. J. 2019, 111, 49–58. [Google Scholar] [CrossRef]
  91. Zhan, M.; Liska, A.J.; Nguy-Robertson, A.L.; Suyker, A.E.; Pelton, M.P.; Yang, H. Data from: Modeled and measured ecosystem respiration in maize–soybean systems over 10 years. Dryad Digit. Repos. 2018. [Google Scholar] [CrossRef]
  92. Gitelson, A. Remote Sensing Estimation of Crop Biophysical Characteristics at Various Scales. In Hyperspectral Remote Sensing of Vegetation; Thenkabail, P.S., Lyon, J.G., Eds.; CRC Press: Boca Raton, FL, USA, 2016; pp. 329–358. ISBN 9780429192180. [Google Scholar]
  93. Gitelson, A.A.; Viña, A.; Verma, S.B.; Rundquist, D.C.; Arkebauer, T.J.; Keydan, G.; Leavitt, B.; Ciganda, V.; Burba, G.G.; Suyker, A.E. Relationship between gross primary production and chlorophyll content in crops: Implications for the synoptic monitoring of vegetation productivity. J. Geophys. Res. 2006, 111, D08S11. [Google Scholar] [CrossRef]
  94. Gitelson, A.A.; Gamon, J.A. The need for a common basis for defining light-use efficiency: Implications for productivity estimation. Remote Sens. Environ. 2015, 156, 196–201. [Google Scholar] [CrossRef] [Green Version]
  95. Wu, L.; Liu, X.; Qin, Q.; Zhao, B.; Ma, Y.; Liu, M.; Jiang, T. Scaling Correction of Remotely Sensed Leaf Area Index for Farmland Landscape Pattern With Multitype Spatial Heterogeneities Using Fractal Dimension and Contextural Parameters. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2018, 11, 1472–1481. [Google Scholar] [CrossRef]
  96. Jin, Z.; Tian, Q.; Chen, J.M.; Chen, M. Spatial scaling between leaf area index maps of different resolutions. J. Environ. Manag. 2007, 85, 628–637. [Google Scholar] [CrossRef]
  97. Ozdogan, M. The spatial distribution of crop types from MODIS data: Temporal unmixing using Independent Component Analysis. Remote Sens. Environ. 2010, 114, 1190–1204. [Google Scholar] [CrossRef]
  98. Gao, F.; Masek, J.; Schwaller, M.; Hall, F. On the blending of the Landsat and MODIS surface reflectance: Predicting daily Landsat surface reflectance. IEEE Trans. Geosci. Remote Sens. 2006, 44, 2207–2218. [Google Scholar] [CrossRef]
  99. Neftel, A.; Spirig, C.; Ammann, C. Application and test of a simple tool for operational footprint evaluations. Environ. Pollut. 2008, 152, 644–652. [Google Scholar] [CrossRef]
  100. Kersebaum, K.C.; Boote, K.J.; Jorgenson, J.S.; Nendel, C.; Bindi, M.; Frühauf, C.; Gaiser, T.; Hoogenboom, G.; Kollas, C.; Olesen, J.E.; et al. Analysis and classification of data sets for calibration and validation of agro-ecosystem models. Environ. Model. Softw. 2015, 72, 402–417. [Google Scholar] [CrossRef]
  101. Bégué, A.; Arvor, D.; Bellon, B.; Betbeder, J.; de Abelleyra, D.; Ferraz, R.P.D.; Lebourgeois, V.; Lelong, C.; Simões, M.R.; Verón, S.; et al. Remote Sensing and Cropping Practices: A Review. Remote Sens. 2018, 10, 99. [Google Scholar] [CrossRef]
  102. Löw, F.; Duveiller, G. Defining the spatial resolution requirements for crop identification using optical remote sensing. Remote Sens. 2014, 6, 9034–9063. [Google Scholar] [CrossRef]
  103. Nguy-Robertson, A.; Gitelson, A.; Peng, Y.; Walter-Shea, E.; Leavitt, B.; Arkebauer, T. Continuous monitoring of crop reflectance, vegetation fraction, and identification of developmental stages using a four band radiometer. Agron. J. 2013, 105, 1769–1779. [Google Scholar] [CrossRef]
  104. Zheng, Y.; Wu, B.; Zhang, M.; Zeng, H.; Zheng, Y.; Wu, B.; Zhang, M.; Zeng, H. Crop Phenology Detection Using High Spatio-Temporal Resolution Data Fused from SPOT5 and MODIS Products. Sensors 2016, 16, 2099. [Google Scholar] [CrossRef]
  105. Onojeghuo, A.O.; Blackburn, G.A.; Wang, Q.; Atkinson, P.M.; Kindred, D.; Miao, Y. Rice crop phenology mapping at high spatial and temporal resolution using downscaled MODIS time-series. GISci. Remote Sens. 2018, 55, 659–677. [Google Scholar] [CrossRef] [Green Version]
  106. Thomason, W.E.; Phillips, S.B.; Davis, P.H.; Warren, J.G.; Alley, M.M.; Reiter, M.S. Variable nitrogen rate determination from plant spectral reflectance in soft red winter wheat. Precis. Agric. 2011, 12, 666–681. [Google Scholar] [CrossRef]
  107. Barnes, E.M.; Clarke, T.R.; Richards, S.E.; Colaizzi, P.D.; Haberland, J.; Kostrzewski, M.; Waller, P.; Choi, C.; Riley, E.; Thompson, T.; et al. Coincident Detection of Crop Water Stress, Nitrogen Status and Canopy Density Using Ground-Based Multispectral Data. In Proceedings of the Fifth International Conference on Precision Agriculture; Robert, P.C., Rust, R.H., Larson, W.E., Eds.; American Society of Agronomy, Crop Science Society of America, Soil Science Society of America: Madison, WI, USA, 2000; pp. 1–15. [Google Scholar]
  108. Franke, J.; Menz, G. Multi-temporal wheat disease detection by multi-spectral remote sensing. Precis. Agric. 2007, 8, 161–172. [Google Scholar] [CrossRef]
  109. Holman, F.; Riche, A.; Michalski, A.; Castle, M.; Wooster, M.; Hawkesford, M.; Holman, F.H.; Riche, A.B.; Michalski, A.; Castle, M.; et al. High Throughput Field Phenotyping of Wheat Plant Height and Growth Rate in Field Plot Trials Using UAV Based Remote Sensing. Remote Sens. 2016, 8, 1031. [Google Scholar] [CrossRef]
  110. Mulla, D.J. Twenty five years of remote sensing in precision agriculture: Key advances and remaining knowledge gaps. Biosyst. Eng. 2013, 114, 358–371. [Google Scholar] [CrossRef]
  111. Zhang, C.; Walters, D.; Kovacs, J.M. Applications of Low Altitude Remote Sensing in Agriculture upon Farmers’ Requests—A Case Study in Northeastern Ontario, Canada. PLoS ONE 2014, 9, e112894. [Google Scholar] [CrossRef]
  112. Neumann, K.; Verburg, P.H.; Stehfest, E.; Müller, C. The yield gap of global grain production: A spatial analysis. Agric. Syst. 2010, 103, 316–326. [Google Scholar] [CrossRef]
  113. Kasampalis, D.; Alexandridis, T.; Deva, C.; Challinor, A.; Moshou, D.; Zalidis, G. Contribution of Remote Sensing on Crop Models: A Review. J. Imaging 2018, 4, 52. [Google Scholar] [CrossRef]
Figure 1. Comparison of actual versus Hybrid-Maize modeled (a) LAI and (b) LUECanopy. The color bars represent the number of points at each marker on the scatter plot.
Figure 1. Comparison of actual versus Hybrid-Maize modeled (a) LAI and (b) LUECanopy. The color bars represent the number of points at each marker on the scatter plot.
Remotesensing 11 01928 g001
Figure 2. (a) Comparison of actual versus Hybrid-Maize modeled LUECanopy. The color bar represents the number of points at each marker on the scatter plot. (b) R2 of actual versus Hybrid-Maize modeled LUECanopy and ΔLUECanopy at different levels of smoothing and values of Δ. N = 2384.
Figure 2. (a) Comparison of actual versus Hybrid-Maize modeled LUECanopy. The color bar represents the number of points at each marker on the scatter plot. (b) R2 of actual versus Hybrid-Maize modeled LUECanopy and ΔLUECanopy at different levels of smoothing and values of Δ. N = 2384.
Remotesensing 11 01928 g002
Figure 3. Comparison of retrieved versus actual (a) LAI from LAIGROUND dataset, (b) LAI from FLUX dataset, and (c) LUECanopy from FLUX dataset from LANDSAT measurements via leave-one-site-out cross validation. The color bars represent the number of points at each marker on the scatter plot.
Figure 3. Comparison of retrieved versus actual (a) LAI from LAIGROUND dataset, (b) LAI from FLUX dataset, and (c) LUECanopy from FLUX dataset from LANDSAT measurements via leave-one-site-out cross validation. The color bars represent the number of points at each marker on the scatter plot.
Remotesensing 11 01928 g003aRemotesensing 11 01928 g003b
Figure 4. RMSE of modeled LAI with respect to ground-truth LAI while varying planting date, seed GDD to maturity, and planting density. Leftmost column represents offset from actual seed GDD to maturity in °C in simulation variant, while second leftmost column represents offset from actual planting density in plants/m2 in simulation variant. Header represents offset from actual planting day in days in simulation variant. Color bar at right and color in main panel represents LAI RMSE for each simulation variant determined by column and header. N = 146.
Figure 4. RMSE of modeled LAI with respect to ground-truth LAI while varying planting date, seed GDD to maturity, and planting density. Leftmost column represents offset from actual seed GDD to maturity in °C in simulation variant, while second leftmost column represents offset from actual planting density in plants/m2 in simulation variant. Header represents offset from actual planting day in days in simulation variant. Color bar at right and color in main panel represents LAI RMSE for each simulation variant determined by column and header. N = 146.
Remotesensing 11 01928 g004
Figure 5. RMSE of modeled LAI with respect to MODIS-retrieved LAI while varying planting date, seed GDD to maturity, and planting density, and MODIS EVI2/LAI regression coefficients. Leftmost column represents offset from actual seed GDD to maturity in °C in simulation variant, while second leftmost column represents offset from actual planting density in plants/m2 in simulation variant. Topmost header represents offset from actual planting day in days in simulation variant. Second topmost header represents slope of EVI2/LAI regression coefficients. Third topmost header represents intercept of EVI2/LAI regression coefficients. Color bar at right and color in main panel represents LAI RMSE for each simulation variant determined by column and header. N = 3280.
Figure 5. RMSE of modeled LAI with respect to MODIS-retrieved LAI while varying planting date, seed GDD to maturity, and planting density, and MODIS EVI2/LAI regression coefficients. Leftmost column represents offset from actual seed GDD to maturity in °C in simulation variant, while second leftmost column represents offset from actual planting density in plants/m2 in simulation variant. Topmost header represents offset from actual planting day in days in simulation variant. Second topmost header represents slope of EVI2/LAI regression coefficients. Third topmost header represents intercept of EVI2/LAI regression coefficients. Color bar at right and color in main panel represents LAI RMSE for each simulation variant determined by column and header. N = 3280.
Remotesensing 11 01928 g005
Figure 6. (ac) Retrieved versus HM Modeled ΔLAI at Δ = 3, 6, 9; (d) Retrieved versus HM Modeled LAI. The colorbars represents the number of points at each marker on the scatter plot.
Figure 6. (ac) Retrieved versus HM Modeled ΔLAI at Δ = 3, 6, 9; (d) Retrieved versus HM Modeled LAI. The colorbars represents the number of points at each marker on the scatter plot.
Remotesensing 11 01928 g006
Figure 7. |r| correlation of the MODIS EVI2 measurements versus HM modeled ΔLUECanopy at different levels of smoothing and values of Δ. N = 2359.
Figure 7. |r| correlation of the MODIS EVI2 measurements versus HM modeled ΔLUECanopy at different levels of smoothing and values of Δ. N = 2359.
Remotesensing 11 01928 g007
Figure 8. Actual versus theoretical |r| performance of the ΔLUECanopy retrievals at all sites at different levels of smoothing and values of Δ. N = 5071.
Figure 8. Actual versus theoretical |r| performance of the ΔLUECanopy retrievals at all sites at different levels of smoothing and values of Δ. N = 5071.
Remotesensing 11 01928 g008
Figure 9. Actual versus theoretical RMSE performance of the ΔLUECanopy retrievals at all sites at different levels of smoothing and values of Δ. N = 5071.
Figure 9. Actual versus theoretical RMSE performance of the ΔLUECanopy retrievals at all sites at different levels of smoothing and values of Δ. N = 5071.
Remotesensing 11 01928 g009
Table 1. Examples of common G × E × M factors included in crop growth model simulations [2,3].
Table 1. Examples of common G × E × M factors included in crop growth model simulations [2,3].
Genotype (G)Environment (E)Management (M)
-Relative maturity/Growing degree days (GDD) to maturity
-GDD to flowering
-Potential kernel number per ear
-Grain growth rate
-Air temperature
-Precipitation
-Solar radiation
-Soil bulk density
-Soil available water
-Soil organic matter
-Soil pH
-Planting date
-Planting density
-Fertilization
-Irrigation
Table 2. Ground-truth data sources.
Table 2. Ground-truth data sources.
NameSource(s)SitesVariables
NameLatitudeLongitudeNameYears
Flux Tower Data (Dataset FLUX)Ameriflux [76]US-Ne1 [35]41.17−96.48GPP
SRAD
Ground-truth LAI
Planting Date
Harvest Date
2001–2009
US-Ne2 [35]41.16−96.472001–2009, odd years
US-Ne3 [35]41.18−96.442001–2009, odd years
US-Ro1 [77]44.71−93.092005, 2009, 2011, 2013
US-Bi2 [78]38.11−121.542017–2018
US-ARM [79]36.61−97.492008
GHG EuropeDE-Kli [80]50.8913.522007, 2012
FR-Gri [81]48.841.952008, 2011
FR-Lam [82]43.51.242006, 2008, 2010
IT-BCi [83]40.5214.962004–2009
NL-Lan [84]51.954.902005
LAI Validation Data (Dataset LAIGROUND)[27]Beltsville39.02−76.85Ground-truth LAI1998
(N = 26)
CEFLES2 [85]44.37–44.460.19–0.412007
(N = 26)
California [86]35.48–39.22−122.14–−119.282011–2012
(N = 59)
Italy (IT-BCi) [83]40.5214.962008–2009
(N = 35)
Mead (US-Ne1 to US-Ne3) [35]41.16−96.462001–2012
(N = 92)
Missouri [87]39.22−92.122002
(N = 10)
NAFE06 [88]−35.08–−34.65145.87–146.32006
(N = 14)
SEN3EXP2009 [85]39.02–39.08−2.13–-2.082009
(N = 10)
SMEX02-IA [89]41.76–42.67−93.73–−93.282002
(N = 21)
SPARC [85]39.03–39.15−2.18–−1.882003–2004
(N = 45)
Table 3. Leave-one-site-out LAIGROUND LANDSAT regression retrieval performance using Equation (4). a and b are the leave-one-site-out regression coefficients defined in Equation (4).
Table 3. Leave-one-site-out LAIGROUND LANDSAT regression retrieval performance using Equation (4). a and b are the leave-one-site-out regression coefficients defined in Equation (4).
Best-Fit CoefficientsLower Bound
Confidence Interval
Upper Bound
Confidence Interval
Site NameLAI RMSENababab
Beltsville0.85268.41−0.927.73−1.188.94−0.65
CEFLES20.60268.55−1.047.76−1.319.10−0.79
California1.32598.19−17.60−1.439.22−0.77
Italy1.58358.49−1.207.82−1.499.33−0.92
Mead1.03927.27−0.715.86−0.97.67−0.03
Missouri0.98108.13−0.877.57−1.188.81−0.64
NAFE060.31148.08−0.857.50−1.429.19−0.61
SEN3EXP20090.89108.20−0.947.61−1.268.90−0.77
SMEX02-IA1.23218.66−1.068.03−1.359.27−0.83
SPARC1.74459.17−1.318.67−1.559.73−1.03
Table 4. Leave-one-site-out FLUX LANDSAT regression retrieval performance using Equations (4) and (5). a, b, c, and d are the leave-one-site-out regression coefficients defined in Equations (4) and (5).
Table 4. Leave-one-site-out FLUX LANDSAT regression retrieval performance using Equations (4) and (5). a, b, c, and d are the leave-one-site-out regression coefficients defined in Equations (4) and (5).
RMSE Best-Fit CoefficientsLower Bound
Confidence Interval
Upper Bound
Confidence Interval
SiteLAILUENabcdabcdabcd
DE-Kli0.850.2049.52−1.241.67−0.169.29−1.361.57−0.209.85−1.111.75−0.13
FR-Gri2.830.1819.52−1.241.67−0.169.28−1.361.58−0.209.88−1.091.76−0.14
FR-Lam1.110.20169.64−1.251.68−0.179.40−1.381.61−0.219.96−1.151.77−0.15
IT-Bci1.410.18329.50−1.271.69−0.179.28−1.391.62−0.229.83−1.151.80−0.15
US-Arm0.140.2319.52−1.241.66−0.169.24−1.361.57−0.199.87−1.031.74−0.13
US-Bi1.630.26129.52−1.251.66−0.169.35−1.401.57−0.209.90−1.171.74−0.13
US-Ne0.830.161248.84−0.801.44−0.095.08−0.961.11−0.189.621.361.680.07
US-Ro1.160.13279.59−1.201.65−0.169.25−1.371.51−0.189.93−1.031.71−0.10
Table 5. Comparison of HM modeled versus retrieved, actual and theoretical |r| and RMSE for the retrieval of ΔLAI and LAI at different values of Δ.
Table 5. Comparison of HM modeled versus retrieved, actual and theoretical |r| and RMSE for the retrieval of ΔLAI and LAI at different values of Δ.
Δ (Days)|r|-Modeled v Retrieved|r|-Modeled v Retrieved
Theoretical
RMSE-Modeled v RetrievedRMSE-Modeled v Retrieved
Theoretical
N
20.520.130.171.462429
30.640.250.291.462429
40.700.360.401.462429
50.750.450.501.462429
60.780.530.591.462429
70.810.590.681.462429
80.830.650.761.462429
90.850.690.841.462429
100.870.730.911.462429
Value Itself (no delta)0.920.881.041.032429
Table 6. Comparison of LAI retrieval performance on all sites except Mead, Nebraska in LAIGROUND dataset trained with actual and HM-modeled Mead, Nebraska LAI values. Only sites with ≥10 points listed site-by-site; all points included in last row.
Table 6. Comparison of LAI retrieval performance on all sites except Mead, Nebraska in LAIGROUND dataset trained with actual and HM-modeled Mead, Nebraska LAI values. Only sites with ≥10 points listed site-by-site; all points included in last row.
SiteNRMSE Trained with Actual DataRMSE Trained with Modeled Data
Beltsville260.840.97
CEFLES2260.770.87
California591.401.39
Italy351.391.26
Missouri100.620.78
NAFE06140.510.47
SEN3EXP2009100.870.79
SMEX02-IA211.201.32
SPARC451.871.83
All except Mead, Nebraska2671.301.29
Table 7. Comparison of LUECanopy retrieval performance on all sites except Mead, Nebraska in FLUX dataset trained with actual and HM-modeled Mead, Nebraska LUECanopy values.
Table 7. Comparison of LUECanopy retrieval performance on all sites except Mead, Nebraska in FLUX dataset trained with actual and HM-modeled Mead, Nebraska LUECanopy values.
SiteNRMSE Trained with Actual DataRMSE Trained with Modeled Data
DE-Kli40.200.20
FR-Gri10.200.10
FR-Lam160.210.29
IT-BCi320.190.35
US-ARM10.220.37
US-Bi2120.260.30
US-Ro1270.130.28
All except Mead, Nebraska930.190.31

Share and Cite

MDPI and ACS Style

Levitan, N.; Kang, Y.; Özdoğan, M.; Magliulo, V.; Castillo, P.; Moshary, F.; Gross, B. Evaluation of the Uncertainty in Satellite-Based Crop State Variable Retrievals Due to Site and Growth Stage Specific Factors and Their Potential in Coupling with Crop Growth Models. Remote Sens. 2019, 11, 1928. https://doi.org/10.3390/rs11161928

AMA Style

Levitan N, Kang Y, Özdoğan M, Magliulo V, Castillo P, Moshary F, Gross B. Evaluation of the Uncertainty in Satellite-Based Crop State Variable Retrievals Due to Site and Growth Stage Specific Factors and Their Potential in Coupling with Crop Growth Models. Remote Sensing. 2019; 11(16):1928. https://doi.org/10.3390/rs11161928

Chicago/Turabian Style

Levitan, Nathaniel, Yanghui Kang, Mutlu Özdoğan, Vincenzo Magliulo, Paulo Castillo, Fred Moshary, and Barry Gross. 2019. "Evaluation of the Uncertainty in Satellite-Based Crop State Variable Retrievals Due to Site and Growth Stage Specific Factors and Their Potential in Coupling with Crop Growth Models" Remote Sensing 11, no. 16: 1928. https://doi.org/10.3390/rs11161928

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop